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spacecraft/payload  mass  properties.  The  second  option  is  a  control 
system  that  is  robust  over  a  smaller  range  of  the  possible  mass 
properties  but  capable  of  estimating  the  mass  properties  and 
adapting  to  these  changing  estimates. 

As  a  simple  example  from  another  type  of  control  problem 
with  similar  demands,  consider  a  robot  arm  capable  of  carrying 
objects  weighing  anywhere  from  1  to  1000  lbs.  The  first  approach  to 
designing  a  controller  for  this  robot  arm  may  be  a  design  robust  to 
this  full  range.  Though  such  a  design  may  be  possible,  performance 
of  this  design  may  not  be  as  accurate  or  timely  as  desired  when 
attempting  to  control  the  arm  with  an  object  near  one  of  the  weight 
extremes.  A  more  attractive  approach  may  be  a  design  robust  only 
to  a  range  of  about  20  lbs  but  capable  of  estimating  the  weight  to 
well  within  a  few  pounds.  Suppose  this  approach  were  used  for 
control  of  a  robot  arm  carrying  a  mass  of  953  lbs.  If  the  control 
system  could  estimate  the  weight  at  950  lbs,  it  cx)uld  adapt  to  this 
estimate  and  provide  robust  control  over  a  20  lb  range  centered  at 
this  estimate.  With  the  control  system  required  to  be  robust  over  a 
much  smaller  range,  it  effectively  has  a  more  precise  model  of  the 
process  to  be  controlled  and  can,  therefore,  provide  more  precise  and 
timely  control. 

This  thesis  employs  the  strategy  of  combining  parameter 
identification  with  robustness  in  design  of  an  attitude  control  system 
for  a  spacecraft  with  a  wide  range  of  mass  properties  performing  an 
aerocapture  maneuver  in  an  environment  with  many  unknowns. 
The  control  law,  the  first  element  of  the  attitude  control  system,  is 
robust  to  a  range  of  system  unknowns  smaller  than  the  full  possible 


Chapter  1 


Introduction 

1.1  Motivation 

When  designing  a  control  system  for  a  particular  process,  a 
control  system  engineer  must  work  with  a  mathematical  model  of 
that  process.  If  the  model  cannot  accurately  predict  the  behavior  of 
the  process,  design  of  an  adequate  control  system  may  not  be 
possible.  The  term  "robustness”  refers  to  a  control  system's 
exhibited  degree  of  immunity  to  system  uncertainties  or  changes. 
Control  systems  can  be  designed  with  a  degree  of  robustness  to  these 
system  uncertainties.  A  tradeoff,  however,  usually  exists  between 
the  amount  of  robustness  and  the  performance  of  the  control  system. 
Control  systems  designed  to  be  robust  over  a  great  range  of  system 
uncertainties  typically  can’t  provide  the  desired  accuracy  or 
responsiveness.  Much  more  accurate  and  efficient  control  can  be 
provided  when  there  is  a  small  range  of  system  unknowns. 

Many  times  it  is  not  possible  to  have  a  highly  accurate 
mathematical  model  of  the  system  to  be  controlled.  For  example,  the 
mathematical  model  of  a  spacecraft  capable  of  carrying  varying 
payloads  depends  greatly  on  the  mass  properties  of  the  unknown 
payload.  There  are  two  options  for  design  of  a  control  system  for 
such  a  spacecraft  with  a  wide  range  of  possible  payloads.  The  first  is 
a  control  system  robust  over  the  full  range  of  the  combined 


range.  A  parameter  identification  scheme,  however,  is  employed  to 
determine  changes  in  these  parameters  real-time,  and  the  control 
system  is  capable  of  adapting  to  the  changing  parameter  estimates. 
This  combination  of  robustness  and  parameter  identification 
effectively  extends  the  range  of  system  unknowns  for  which  attitude 
control  is  possible. 

1.2  Mission  Description 

As  a  spacecraft  approaches  a  planet  on  a  hyperbolic  trajectory, 
it  must  effect  a  large  velocity  reduction  (AV)  in  order  to  enter  a 
desired  orbit  about  the  planet.  The  usual  method  for  performing 
such  a  AV  has  been  to  use  rocket  firing.  Unfortunately,  the  required 
firing  forces  the  spacecraft  to  carry  a  large  mass  of  propellent  for 
this  purpose. 

Studies  [1]  have  investigated  using  an  aeroassisted  maneuver, 
called  aerocapture,  to  generate  most  of  the  AV  required  to  obtain  a 
desired  orbit.  With  aerocapture,  the  spacecraft  descends  into  the 
upper  portions  of  the  planetary  atmosphere  as  it  passes  by  and  uses 
the  aerodynamic  forces  generated  to  reduce  velocity.  Ideally,  the 
vehicle  exits  the  atmosphere  with  the  velocity  required  to  reach  a 
desired  orbital  altitude.  Once  the  target  altitude  is  reached,  a 
comparatively  small  rocket  firing  is  required  to  obtain  the  desired 
orbit  characteristics  (e.g.,  circularization). 

The  smaller  propulsive  AV  required  for  a  spacecraft  employing 
aerocapture  reduces  the  amount  of  fuel  needed,  allowing  for  a  larger 
payload  or  smaller  vehicle  for  the  same  mission.  Even  considering 
the  added  weight  of  an  aeroshell  and  the  thermal  protection  system 


needed  during  atmospheric  flight,  the  reduced  weight  of  fuel  can 
result  in  as  much  as  twice  the  payload  capacity  of  a  similar  all- 
propulsive  vehicle  [2].  This  benefit  has  made  aerocapture  a  likely 
element  of  many  future  missions. 

Aerocapture  at  Mars  has  been  advocated  for  future  Martian 
missions  along  with  Earth  aerocapture  for  missions  returning  to  Earth 
from  high  Earth  orbits,  the  moon,  and  Mars  [3].  To  underscore  the 
interest  in  aerocapture,  the  December  1990  report  from  the  Advisory 
Committee  on  the  Future  of  the  U.S.  Space  Program  headed  by  Mr. 
Norman  Augustine  recommends  pursuing  the  technology  of 
aerocapture  [4].  Also,  on  March  12,  1991,  Japan's  Hiten  satellite 
became  the  first  spacecraft  to  demonstrate  an  Earth  aerocapture 
after  circling  the  Moon  in  1990  [5]. 

1.3  Attitude  Control  System  Requirements 

The  main  disadvantage  of  an  aerocapture  maneuver  relative  to 
an  all-propulsive  maneuver  is  the  difficulty  in  properly  guiding  and 
controlling  the  vehicle  while  in  the  atmosphere.  Errors  in  the  final 
orbit  are  highly  sensitive  to  atmospheric  exit  conditions.  Small 
errors  in  the  exit  velocity  or  flight  path  angle  could  cause  large 
errors  in  the  final  orbit.  Extra  fuel  would  then  be  needed  to  correct 
to  the  desired  orbit,  reducing  the  original  benefit  of  the  aerocapture 
maneuver.  Aerocapture  guidance  and  control  systems  must  be 
robust  to  a  wide  range  of  uncertainties  to  avoid  these  exit  errors. 

In  simple  terms,  guidance  determines  where  a  vehicle  must  go 
in  order  to  meet  mission  requirements,  while  control  determines 
what  the  vehicle  must  do  to  get  there.  For  the  purposes  of  this 


thesis,  guidance  is  considered  to  be  the  determination  of  the 
commands  needed  to  follow  the  trajectory  required  to  meet  the 
objectives  of  the  aerocapture  mission.  Since  the  trajectory  is 
dependent  on  the  aerodynamic  lift  and  drag  acting  on  the  vehicle, 
guidance  consists  of  the  calculation  of  the  vehicle  attitude  needed  to 
generate  the  required  aerodynamic  forces.  Attitude  control  is  the 
determination  of  how  to  employ  the  spacecraft's  control  actuators 
(e.g.,  jets,  aerosurfaces,  and  control  moment  gyros)  in  a  manner  to 
achieve  the  attitude  commanded  by  guidance.  This  thesis  addresses 
the  attitude  control  problem. 

An  aerocapture  guidance  law  must  be  capable  of  gui^^lng  a 
vehicle  through  an  atmosphere  to  ensure  it  exits  the  atmosphere 
with  the  proper  energy  to  reach  a  desired  orbit.  Uncertainties  that  a 
guidance  law  must  accommodate  include  atmospheric  density 
fluctuations  and  aerodynamic  prediction  uncertainties.  At  Earth, 
density  fluctuations  of  up  to  30%  off  the  nominal  values  can  be 
expected  [6].  At  planets  with  lesser  known  atmospheres,  such  as 
Mars,  larger  fluctuations  are  considered  possible.  Little  data  is 
currently  available,  but  densities  ranging  from  -50%  to  -(-100%  off 
present  Martian  density  models  are  considered  likely.  Faulty 
guidance  could  cause  a  spacecraft  to  exit  the  atmosphere  without  the 
proper  energy  to  reach  the  target  orbit  altitude.  At  one  extreme, 
diving  too  deep  could  reduce  the  energy  to  the  point  that  the 
spacecraft  would  not  be  able  to  exit,  sending  the  vehicle  crashing  into 
the  planet’s  surface.  At  the  other  extreme,  going  too  shallow  might 
cause  the  spacecraft  to  skip  out  of  the  atmosphere  with  too  much 
energy  to  be  captured  into  orbit  about  the  planet. 


Most  current  aerocapture  guidance  algorithms  use  vehicle  bank 
angle  0  as  the  only  control  (angle-of-attack,  a,  and  sideslip  angle,  (3, 
are  assumed  constant)  [2,7,8,91.  For  the  purposes  of  this  thesis,  this 
set  of  angles  (i.e.,  0,  a,  and  P)  is  called  the  "velocity  angles".  These 
angles  are  used  to  describe  the  attitude  of  the  vehicle  with  respect  to 
a  frame  containing  the  velocity  vector  as  defined  in  Appendix  A.  In¬ 
plane  guidance  generates  bank  angle  commands  that  vary  the 
direction  of  the  lift  vector  in  order  to  fly  a  reference  atmospheric 
trajectory  which  allows  the  spacecraft  to  reach  the  desired  orbital 
altitude.  Out-of-plane  guidance  periodically  commands  bank  angle 
reversals  in  order  to  use  out-of-plane  lift  to  maintain  the  vehicle  in 
the  desired  orbital  plane.  Future  guidance  systems  may  also  vary 
angle-of-attack  to  regulate  the  total  amount  of  lift,  allowing  in-plane 
guidance  to  be  performed  without  generating  any  undesired  out-of¬ 
plane  lift.  For  the  purpose  of  this  thesis,  only  bank  modulation  will 
be  considered. 

The  task  of  the  attitude  control  system,  then,  is  to  implement 
the  bank  commands  and  periodic  reversals  while  maintaining  trim  a 
and  p.  The  attitude  control  system  must  perform  this  task  in  the 
presence  of  uncertainties.  Variations  in  the  vehicle  mass  properties 
(e.g.,  total  mass,  inertia,  and  location  of  the  center  of  mass)  can  have 
a  great  impact  on  the  effectiveness  of  control  actuators.  For  generic 
"bus"  vehicles  capable  of  carrying  and  delivering  several  types  of 

payloads,  control  must  be  performed  over  a  wide  range  of  mass 

properties.  Rather  than  develop  or  modify  a  control  system  for  each 

payload,  a  design  capable  of  controlling  over  this  range  of  mass 

properties  is  highly  desirable.  In  addition,  properties  of  the  control 
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actuators  themselves  (e.g.,  jet  thrust  levels  and  directions, 
aerodynamic  properties  of  control  aerosurfaces)  and  atmospheric 
conditions  (e.g.,  density)  can  vary,  adding  to  the  uncertainties  with 
which  the  control  system  must  contend. 

1.4  Prior  Work 

Several  guidance  laws  have  been  developed  for  two  particular 
aerocapture  missions  [2,7, 8,9].  Considerable  work  has  gone  into 
developing  a  guidance  system  for  NASA’s  Aeroassist  Flight 
Experiment  (AFE).  The  AFE,  expected  to  be  flown  in  the  later  1990's, 
consists  of  a  small,  blunt  (lift-to-drag  ratio  =  0.3)  vehicle,  sketched  in 


Figure  1.1.  AFE  Vehicle 
2  1 


Figure  1.1,  to  be  deployed  from  the  Space  Shuttle  with  the  objective 
of  simulating  a  vehicle  transferring  from  geosynchronous  Earth  orbit 
to  a  low  Earth  orbit  using  aerocapture.  Once  deployed  from  the 
Shuttle's  cargo  bay,  the  AFE  vehicle  will  be  propelled  into  the  Earth’s 
atmosphere  by  a  solid  rocket  motor  (SRM)  burn,  achieving  entry 
interface  conditions  matching  those  of  a  vehicle  returning  from 
geosynchronous  Earth  orbit.  After  the  atmospheric  pass,  the  AFE 
orbit  will  be  circularized  with  another  rocket  burn,  allowing  the 
vehicle  to  be  retrieved  by  the  Space  Shuttle  [10].  Figure  1.2  outlines 
the  main  events  during  the  AFE  flight. 

The  other  mission  that  has  received  considerable  attention  is 
the  Mars  Rover  Sample  Return  (MRSR)  mission,  studied  during  1986- 
89  [11,12,13,14,15].  In  this  mission,  an  orbiter,  lander,  and  rover, 
contained  in  an  aeroshell,  would  aerocapture  into  a  Martian  orbit 
from  which  a  landing,  sample  collection,  and  sample  return  to  Earth 
would  be  accomplished.  Since  the  president's  speech  on  July  20, 
1989  calling  for  a  human  expedition  to  Mars,  the  MRSR  mission  has 
been  set  aside.  It  is  likely,  however,  that  any  resulting  Martian 
exploration  missions  will  include  many  of  the  elements  of  the  MRSR 
mission,  including  aerocapture  [4]. 

While  there  has  been  considerable  effort  in  the  development  of 
aerocapture  guidance  laws,  little  effort  has  been  focused  on 
aerocapture  attitude  control  system  design.  A  simple  system 
consisting  of  a  proportional  control  law  and  a  table  look-up  Jet 
selection  procedure  has  been  developed  for  the  AFE  mission  [16]. 
This  attitude  control  system  is  adequate  for  the  one-time,  tightly 
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constrained  AFE  experiment  where  vehicle  mass  and  actuator 
properties,  as  well  as  initial  entry  conditions,  are  well  defined.  Such 
a  simple  attitude  control  system,  however,  cannot  provide  the 
robustness  required  for  a  mission  with  vehicle  and  environmental 
properties  varying  over  a  wide  range. 


Figure  1.2.  AFE  Mission 
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1.5  Thesis  Goal 


1.5.1  Goal 

The  goal  of  this  thesis  is  an  attitude  control  system  design 
applicable  to  aerocapture  and  capable  of  meeting  some  of  the 

requirements  outlined  in  Section  1.3  above.  A  desirable  ultimate 
goal  is  a  generic  attitude  control  system  for  use  with  any  aerocapture 
mission  (e.g.,  at  Earth  or  at  Mars,  with  varying  ranges  of  the 

atmospheric  unknowns  affecting  the  vehicle  aerodynamics)  and  any 

vehicle  (e.g.,  high  vs.  low  lift-to-drag  ratio,  with  jets  and/or 
aerosurface  control  actuators,  capable  of  carrying  a  single  payload  or 
a  wide  range  of  payloads).  It  is  also  desired  that  the  system  be  able 
to  control  both  variable  bank  angle  and  angle-of-attack  and  have  the 
capability  to  identify  any  and  all  of  the  following:  failed  jets, 
reduced  or  misaligned  jet  thrust,  changing  vehicle  mass  properties, 
and  changes  in  aerosurface  effectiveness  due  to  the  aerosurface  itself 
or  to  changes  in  atmospheric  conditions. 

This  thesis  addresses  a  subset  of  the  ultimate  goal.  It  develops 
a  candidate  attitude  control  system  capable  of  performing  attitude 
control  during  an  aerocapture  maneuver  for  a  vehicle  with  a  wide 

range  of  possible  mass  properties.  Actuator  properties  (e.g.,  thrust 
levels  and  directions)  are  considered  constant  and  known.  Since  the 
AFE  vehicle  and  mission  are  presently  the  best  defined  for  an 
aerocapture  mission,  the  candidate  control  system  has  been  designed 
for  the  AFE  mission  scenario  (i.e..  Earth  aerocapture)  and  for  a 
vehicle  similar  to  the  actual  AFE  vehicle  but  with  a  wider  range  of 


possible  mass  properties.  It  is  assumed  that  control  is  possible  using 
the  AFE's  reaction  control  system  (RCS)  jets,  the  only  type  of  control 
actuator  employed  for  the  AFE. 


1.5.2  Assumptions 

The  following  is  a  summary  of  the  basic  assumptions  used 
throughout  this  thesis. 

1.  The  control  system  is  designed  to  provide  attitude  control 
for  a  vehicle  performing  an  Earth  aerocapture. 

2.  The  attitude  control  system  must  implement  bank  angle 
commands  while  uiaintaining  trim  angle-of-attack  and 
sideslip. 

3.  Control  to  within  2  degs  of  the  commanded  velocity 
angles  is  desired  [17]. 

4.  The  vehicle  to  be  controlled  is  rigid. 

5.  The  aerocapture  vehicle  to  be  controlled  has  the  same 
control  actuators  (i.e.,  RCS  jets  only)  and  aerodynamic 
properties  as  the  AFE  vehicle. 

6.  The  aerocapture  vehicle  is  capable  of  a  much  wider  range 
of  mass  properties  than  the  AFE  vehicle. 

7.  Jet  locations  and  thrust  vectors  are  constant  and  known. 

8.  Sensors  on  the  vehicle  provide  ideal,  lag-free 
measurements  of  the  vehicle  attitude  and  rates. 

9.  The  aerodynamic  torques  experienced  by  the  vehicle 
during  the  aerocapture  maneuver  are  small  compared  to 
the  torques  applied  by  the  RCS  jets. 
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1.5.3  Intended  Contributions 

The  intended  contributions  of  this  thesis  can  be  summarized  as 
follows: 

1.  The  development  of  a  robust  control  law  for  attitude 
control  of  a  vehicle  with  variable  mass  properties. 

2.  The  combination  of  the  control  law  with  an  adaptable  jet 
selection  algorithm  to  create  a  control  system  robust  to  a 
given  range  of  mass  properties  but  which  is  smaller  than 
the  full  range  possible. 

3.  The  integration  of  a  real-time  mass  property  identiti- 
cation  algorithm  to  extend  the  range  of  mass  properties 
for  which  the  candidate  control  system  is  robust. 

1.6  Thesis  Overview 

The  remaining  chapters  describe  the  candidate  attitude  control 
system  design  and  its  performance  during  computer  simulation 
testing.  Chapter  2  provides  a  system  overview  of  the  candidate 
design.  The  designs  of  the  three  main  elements  of  the  attitude 
control  system,  the  control  law,  the  actuator  selection  algorithm,  and 
the  mass  property  identification  algorithm,  are  outlined  in  Chapters 
3,  4,  and  5,  respectively.  Chapter  6  then  describes  the  integration  of 
the  main  elements  into  a  functional  attitude  control  system. 

Testing  of  the  candidate  design  is  accomplished  using  a  6 
degree-of-freedom  (6-DOF)  FORTRAN  computer  simulation  of  an 
Earth  aerocapture.  This  simulation  and  the  testing  results  are 
discussed  in  Chapter  7. 
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The  conclusions  and  contributions  drawn  from  this  research  are 
summarized  in  Chapter  8  along  with  recommendations  for  future 
research. 
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Chapter  2 

System  Overview 


2.1  Control  System  Structure 

The  candidate  attitude  control  system  consists  of  three  main 

elements:  the  control  law,  the  actuator  selection  algorithm,  and  the 
parameter  identification  algorithm  as  depicted  in  Figure  2.1. 

The  first  element,  the  3-dimensional  attitude  control  law, 
compares  the  commanded  velocity  angles  from  the  guidance  law  to 
the  estimated  velocity  angles  and  determines  the  angular 
accelerations  required  to  implement  the  guidance  commands.  It  is 
designed  to  exhibit  only  the  degree  of  robustness  required  for 
control  over  the  range  of  possible  values  of  the  inertia  matrix  for  the 
actual  AFE  vehicle  as  presented  in  Appendix  B. 

The  required  angular  accelerations  from  the  control  law  are  the 
commanded  inputs  to  the  second  main  element,  the  actuator 

selection  algorithm.  This  algorithm  uses  the  control  system's 

knowledge  of  the  vehicle  properties  to  determine  the  effectiveness  of 
each  actuator,  and  then  calculates  the  actuator  activity  required  to 
generate  the  commanded  accelerations.  Jet  duty  cycles  and,  for  a 
more  general  design,  aerosurface  deflections  are  commanded  by  the 
actuator  selection  and  implemented  by  the  associated  actuators. 

The  vehicle  dynamics  are  determined  by  the  moments  and 
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Figure  2.1.  Aerocapture  Control  System  Structure 
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forces  acting  on  the  vehicle.  These  moments  and  forces  consist  of 
those  applied  by  the  control  system  and  those  from  the  environment. 

Lag-free  measurements  of  the  vehicle's  attitude  and  rates  are 
made  by  sensors  in  the  vehicle’s  inertial  measurement  unit  (IMU). 
The  measured  attitude  is  used  to  generate  the  velocity  angle 
estimates  used  by  the  control  law  while  the  measured  angular  rates 
are  used  by  the  third  main  element  of  the  control  system,  the  mass 
property  identification  algorithm.  This  algorithm  compares  the 
vehicle’s  angular  rate  changes,  determined  from  the  IMU  rate 
measurements,  to  the  anticipated  rate  changes  due  to  actuator 
activity  predicted  by  models  of  the  vehicle  and  its  control  actuators. 
This  comparison  is  used  to  produce  estimates  of  unknown  and 
changing  mass  properties.  Updated  estimates  are  then  supplied  to 
the  control  law  and  the  actuator  selection  algorithm,  increasing  the 
range  of  mass  properties  for  which  the  attitude  control  system  is 
usable. 

The  major  elements  of  the  attitude  control  system  developed 
here  have  been  adapted  from  previous  works.  The  control  law,  the 
first  main  element,  is  a  robust,  nonlinear  design  based  on  the  sliding 
mode  control  techniques  outlined  in  references  [18,19].  A  variation 
of  the  actuator  selection  algorithm  used  previously  for  blended 
actuator  selection  [20,21]  is  used  to  create  a  jet  selection  algorithm 
capable  of  adjusting  to  changing  mass  property  estimates.  Mass 
property  identification  is  performed  using  the  second  order, 
nonlinear  filter  design  developed  in  references  [22,23,24].  The  main 
contribution  of  this  thesis  is  the  adaptation  and  integration  of  the 


existing  elements  to  form  a  viable  attitude  control  system  design 

with  application  to  aerocapture. 

2.2  Control  Law 

The  candidate  control  system  requires  a  control  law  that  can 

determine  the  angular  accelerations  needed  by  the  vehicle  to  track 

the  velocity  angles  commanded  by  the  guidance.  Since  the  vehicle 
dynamics  are  affected  by  variable  vehicle  mass  properties  and 

atmospheric  conditions,  the  control  law  must  be  robust  to  the 

anticipated  ranges  of  these  parameters. 

Linear  control  design  strategies  such  as  Hoo,  H2,  and  \L  synthesis 
[25,26,27]  can  be  used  for  designing  robust  control  laws.  Sliding 
mode  control,  however,  was  chosen  for  this  design.  The  dynamics  of 
a  spacecraft  are  inherently  nonlinear.  With  sliding  mode  control,  an 
extensive,  linearized  system  model  is  not  required.  As  long  as 

bounds  on  parametric  unknowns  and  unmodelled  dynamics  are 

known,  a  relatively  simple,  nonlinear  model  is  all  that's  required. 
Sliding  mode  control  then  provides  precise  control  robust  over  the 

range  of  bounded  unknowns. 

2.3  Actuator  Selection 

Once  the  control  law  calculates  the  angular  accelerations 
required  for  the  vehicle  to  implement  the  guidance  commands,  the 
actuator  selection  algorithm  must  determine  which  control  actuators 
to  activate,  and  to  what  extent,  in  order  to  generate  the  commanded 
accelerations.  With  only  RCS  jets  available  to  the  candidate  system, 
the  required  jet  selection  algorithm  must  determine  which  jets  to 
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fire,  and  with  what  duty  cycle,  to  give  the  aerocapture  vehicle  the 
commanded  angular  accelerations. 

In  order  for  an  actuator  selection  algorithm  to  perform  its  task, 
it  must  have  accurate  knowledge  of  the  effectiveness  of  each  of  its 
control  actuators.  The  effectiveness  of  an  actuator  depends  on  many 
things  including  the  location  of  the  center  of  mass,  the  moments  of 
inertia  of  the  vehicle,  the  location  and  orientation  of  the  actuator,  and 
the  level  of  activity  (e.g.,  amount  of  thrust  of  a  jet)  provided  by  the 
actuator. 

For  the  one-time  AFE  mission  with  well  defined  mass  and 
actuator  properties,  the  effectiveness  of  each  jet  is  calculated  a  priori 
and  jet  selection  is  performed  via  table  look-up.  For  example,  if  the 
commanded  roll  acceleration  is  -»-5  deg/sec^,  then  jets  #1  and  #2  are 
fired  during  the  next  time  step.  If  the  commanded  yaw  acceleration 
is  -0.5  deg/sec2,  then  jet  #8  is  fired  (see  Figure  4.1  in  Chapter  4  for 
jet  numbering  and  locations).  A  table  look-up  selection  algorithm 
like  this  which  uses  a  priori  effectiveness  calculations  cannot  always 
perform  its  task  of  implementing  commanded  angular  accelerations 
when  vehicle  properties  are  not  initially  known  or  change  during  a 
mission.  For  the  AFE  vehicle  with  nominal  moment  arms  of  a  little 
under  4  feet,  a  2  foot  change  in  the  center  of  mass  location  greatly 
changes  the  effectiveness  of  the  RCS  jets.  Without  a  method  of 
identifying  such  a  change,  the  table  look-up  jet  selection  algorithm 
must  choose  jet  firings  based  on  erroneous  measures  of  jet 
effectiveness.  Control  would  be  inefficient  at  best,  with  the 
possibility  of  being  unstable.  For  a  mission  where  vehicle  properties 
are  unknown  a  priori  or  can  change  significantly,  an  actuator 


selection  algorithm  that  can  adapt  to  the  changing  properties  during 
the  mission  is  required. 

Several  methods  are  available  for  design  of  an  adaptable 
actuator  selection  algorithm.  Included  among  those  available  are  the 
dot  product,  control  axes,  pseudo  inverse,  and  linear  programming 
methods  outlined  in  reference  [28].  Being  the  most  fuel  optimal  and 
adaptable  of  these  methods,  the  linear  programming  approach  was 
chosen  to  perform  the  actuator  selection  task  for  this  thesis.  This 
approach  also  has  the  benefit  of  having  already  been  successfully 
flight  tested  in  an  experimental  jet  selection  algorithm  for  the  Space 
Shuttle  on  missions  STS-51G  (June  1985)  and  STS-61B  (November 
1985). 


2.4  Parameter  Identification 

In  order  for  the  actuator  selection  algorithm  to  adapt  to 
changing  vehicle  parameters  and  to  increase  the  overall  robustness 
of  the  attitude  control  system,  a  parameter  identification  algorithm 
capable  of  estimating  and  updating  changing  or  unknown  vehicle 
parameters  is  required. 

For  this  thesis,  the  jet  properties,  such  as  thrust  level  and 
direction,  are  considered  fixed  and  known,  but  the  vehicle  mass 
properties  (i.e.,  the  vehicle  inertia  matrix  and  location  of  the  center 
of  mass)  are  assumed  to  be  initially  unknown  and  subject  to  change. 
Mass  property  identification  and  estimation  is  accomplished  using  a 
second-order,  nonlinear  filter  design  resembling  an  extended  Kalman 
filter  [23,24].  The  algorithm  uses  a  model  of  the  dynamics  of  a  rigid 
spacecraft  to  predict  the  output  of  rate  gyros  and  accelerometers  due 
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to  jet  firings.  Comparing  its  predictions  to  the  measured  output, 

revisions  of  mass  property  estimates  are  made. 

This  algorithm  will  be  used  to  estimate  the  AFE  mass 

properties  using  data  gathered  in  flight.  These  estimates  will  be 

compared  to  actual  measurements  of  the  mass  properties  made  on 
the  ground  prior  to  and  after  the  AFE  mission.  Because  of  this 

association  with  the  AFE,  this  algorithm  is  the  chosen  approach  for 
mass  property  identification  for  this  problem. 
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Chapter  3 


Control  Law  Design 

3.1  Introduction 

The  concepts  of  sliding  mode  control  are  used  in  designing  a 
robust,  nonlinear  control  law  which  determines  the  angular 
accelerations  required  to  track  the  bank  angle  trajectory  commanded 
by  the  guidance  law  while  maintaining  trim  a  and  p.  This  chapter 
covers  the  development  of  this  first  element  of  the  candidate 
attitude  control  system.  The  nonlinear  stability  analysis  theory 
which  provides  the  basis  for  sliding  mode  control  is  presented  in 
Section  3.2.  Section  3.3  develops  the  basic  sliding  mode  concepts 
through  the  design  of  a  control  law  for  a  simple  single-input,  single¬ 
output  system.  A  further  example  of  sliding  mode  control  design  for 
the  single-input,  single-output  case  is  outlined  in  Section  3.4,  while 
Section  3.5  describes  the  design  of  the  sliding  mode  attitude  control 
law  for  the  multi-input,  multi-output  aerocapture  vehicle. 

3.2  Fundamentals  of  Lyapunov  Theory 

Sliding  mode  control  originated  from  the  stability  theory  for 
nonlinear  systems  introduced  by  the  Russian  mathematician  A.M. 
Lyapunov  in  the  late  1800's.  Chapter  3  of  reference  [19]  and  chapter 
3  of  reference  [29]  provide  good  outlines  of  the  fundamentals  of 
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Lyapunov  theory.  The  following  ii>  a  summary  of  the  main  points 
presented  in  these  references. 

3.2.1  Stability 

Before  discussing  Lyapunov  stability  theory,  the  basic 
definitions  of  stability  first  need  to  be  discussed.  Consider  the  basic 
nonlinear  system  of  the  form 

x(t)=f(x,  t)  (3.1) 

An  equilibrium  state  Xeq  of  this  system  is  a  state  where  once  x(t)  = 
Xcq,  it  remains  equal  to  Xeq  for  all  time.  In  mathematical  terms,  an 
equilibrium  state  satisfies  0  =  f(^q,t)  for  all  time.  An  equilibrium 
state  is  considered  stable  if  the  system  trajectory  remains  arbitrarily 
close  to  the  equilibrium  state  when  starting  sufficiently  close  to  it. 
Asymptotic  stability  implies  not  only  that  an  equilibrium  state  is 
stable,  but  also  that  system  trajectories  starting  sufficiently  close  to 
the  equilibrium  state  actually  converge  to  the  equilibrium  state  as 
time  goes  to  infinity.  An  equilibrium  state  which  is  stable  but  not 
asymptotically  stable  is  often  called  marginally  stable.  Finally,  global 
asymptotic  stability  of  an  equilibrium  state  means  that  asymptotic 
stability  holds  for  any  initial  state. 

3.2.2  Lyapunov's  First  Method 

Lyapunov  presented  two  methods  for  analyzing  stability  of 
nonlinear  systems.  With  the  first  method  of  Lyapunov,  the  stability 
of  an  equilibrium  state  of  a  nonlinear  system  of  the  form  of  Equation 
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3.1  can  be  analyzed  by  first  linearizing  the  system  about  the 
equilibrium  state.  Conclusions  as  to  the  stability  of  the  nonlinear 
system  can  be  drawn  based  on  the  stability  of  the  linearized  sys  cm. 
If  the  linearized  system  is  strictly  stable,  the  nonlinear  system  is 
stable  at  the  equilibrium  state.  If  the  linearized  system  is  unstable, 
the  equilibri  m  stat“  is  also  unstable.  If  the  linearized  system  is 
marginally  stable,  no  conclusions  can  be  drawn  as  to  the  stability  of 
the  nonlinear  system  at  the  equilibrium  state  [19]. 

3.2.3  Lyapunov's  Second  Method 

The  basic  approach  of  the  second  method  of  Lyapunov,  also 
called  Lyapunov's  direct  method,  is  to  try  to  find  a  scalar  function, 
based  on  the  nonlinear  system  in  question,  that  meets  the  criteria  of 
a  Lyapunov  function.  To  be  a  Lyapunov  function  of  a  system,  the 
function  must  exhibit  two  properties.  The  first  is  that  the  function 
V(x)  must  be  positive  definite.  For  V(x)  to  be  positive  definite,  the 


following 

must  be  true: 

(a) 

V(x)  must  have  continuous  partial 
respect  to  the  components  of  x 

derivatives 

with 

(b) 

V(0)  =  0 

(c) 

V(x)  >  0  for  X  0 

The  second  property  required  of  a  Lyapunov 

function  is  the 

derivative 

of  V(x)  must  be  negative  definite  inside 

a  region  R 

about 

the  origin. 

To  be  negative  definite,  the  same  conditions  as  those 

for  a 

positive  definite  function  must  apply  but  with  the  inequality  sign  in 
condition  (c)  reversed. 

Given  these  properties  of  a  Lyapunov  function,  the  theorem 
used  in  Lyapunov's  second  method  of  stability  analysis  simply  states 
that  the  null  solution  to  Equation  3.1  (i.e.,  the  origin)  is 
asymptotically  stable  in  a  region  R  around  the  origin  if  there  exists  a 
Lyapunov  function  of  the  system  over  the  region  R  possessing  the 

two  properties  descri'o>ed  above.  Further,  if  the  function  V(x) 
satisfies  the  conditions  of  a  Lyapunov  function  for  all  x,  the 
equilibrium  state  at  the  origin  is  globally  asymptotically  stable.  This 
theorem  provides  only  sufficient  conditions  for  stability 

determination.  Failure  to  find  a  Lyapunov  function  does  not  prove 
instability  but  only  represents  failure  to  prove  stability  [29]. 

The  difficulty  in  applying  Lyapunov  theory  to  stability  analysis 
of  nonlinear  systems  is  knowing  what  to  choose  as  a  possible 
Lyapunov  function.  For  a  given  system,  there  may  be  many  possible 

Lyapunov  functions  of  the  system;  only  one  is  needed  to  determine 

stability.  The  following  example  points  to  a  typical  source  of  v 
Lyapunov  function  for  a  mechanical  system  -  the  total  mechanical 
energy. 

3.2.4  Lyapunov  Stability  Analysis  of  a  Simple  Pendulum 

As  an  illustration  of  the  second  method  of  Lyapunov,  consider 
the  pendulum  in  Figure  3.1  [19].  The  dynamics  of  this  simple 
pendulum  are  given  by  the  nonlinear  equation 

mr^Q  +  kO  +  mgr  (sin  6)  =  0  (3.2) 
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where  m  is  the  mass  of  the  pendulum,  r  is  the  length  of  the 
pendulum,  g  is  the  acceleration  due  to  gravity,  k  is  the  coefficient  of 
friction,  and  0  is  the  angle  defined  in  Figure  3.1. 


Figure  3.1.  Simple  Pendulum 


Defining  the  state  vector  as  x  =  [  0  0  ]T,  a  stable  equilibrium 
state  of  this  system  is  at  Xeq  =  [0  0]^  representing  the  pendulum  at 
rest  at  the  bottom  of  its  arc.  Starting  the  pendulum  at  this  point  with 
no  velocity,  the  pendulum  would  stay  at  this  point.  Intuitively,  this 
point  is  also  asymptotically  stable  for  all  starting  trajectories  other 
than  at  x(0)  =  [jt  0]^'  which  represents  the  inverted  pendulum  at  rest. 
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an  unstable  equilibrium  state.  As  the  pendulum  is  raised  from  0=0, 
its  potential  energy  is  increased.  When  released,  the  pendulum 
swings  back  towards  0=0,  trading  potential  energy  for  kinetic 
energy  until  it  passes  the  bottom  of  its  arc  and  begins  swinging  back 
upward,  trading  kinetic  energy  for  potential  energy.  For  a  pendulum 
with  no  dissipative  friction  (i.e.,  k  =  0),  the  total  mechanical  energy  of 
the  system  (i.e.,  kinetic  +  potential)  remains  constant  with  the 
pendulum  swinging  back  and  forth,  always  returning  to  the  height 
from  which  it  was  released.  With  friction,  however,  the  total 
mechanical  energy  of  the  pendulum  is  dissipated  and  each  upward 
swing  of  the  pendulum  will  not  reach  the  height  of  its  previous 
swing.  Eventually,  the  pendulum  will  slow  down  and  come  to  rest  at 
the  bottom  of  its  arc,  representing  the  equilibrium  point  Xeq  =  [0  0]^. 

Before  coming  to  rest,  the  pendulum  possesses  positive  total 
mechanical  energy.  Since  the  friction  causes  the  total  mechanical 
energy  to  dissipate,  the  total  mechanical  energy  can  be  thought  of  as 
having  a  negative  derivative.  These  two  characteristics  of  the  total 
mechanical  energy  (i.e.,  positive  with  negative  derivative)  are  the 
desired  characteristics  of  a  Lyapunov  function  suggesting  the 
possibility  of  a  Lyapunov  function  based  on  the  total  mechanical 
energy  of  the  system. 

The  total  mechanical  energy  of  the  pendulum  is  given  by  the 


function 


representing  the  combined  kinetic  and  potential  energy.  Applying 
the  second  method  of  Lyapunov,  V(x)  meets  the  two  criteria 
required  of  a  Lyapunov  function  for  all  states  except  x  =  [i:  0]^: 

1.  V(x)  is  positive  definite 

(a)  V(x)  has  continuous  partial  derivatives  with 
respect  to  the  components  of  x 

(b)  V(0)  =  0 

(c)  V(x)  >  0  for  X  56  0,  [ti  0]*^ 

2.  Differentiating  V(x)  and  substituting  in  for  0  from 
Equation  3.2,  V(x)  is  shown  to  be  negative  definite: 

>^x)  =  mr^  00  +  mgr  sin  0  0 

=  mr^  ( — ^  0  -  ^  sin  01 0  +  mgr  sin  0  0 
I  mr2  r  ) 

V(x)  =  -k  0^  <  0,  0  0 

The  function  V(x)  is  a  Lyapunov  function  of  the  simple  pendulum  for 
all  X  except  x  =  [tc  0]^,  For  the  simple  pendulum,  the  equilibrium 
point  at  the  origin  is  asymptotically  stable  everywhere  except  for  the 
single  state  x  =  [n  0]T,  the  result  expected  from  intuition. 

3.3  Fundamentals  of  Sliding  Mode  Control 

This  section  outlines  the  basics  of  sliding  mode  control  as 
presented  in  chapter  7  of  reference  [19].  The  concept  of  robust 
sliding  mode  control  stems  from  the  second  method  of  Lyapunov 
outlined  above.  The  basic  idea  is  to  design  a  feedback  control  law 
such  that  a  function  based  on  the  closed  loop  system  can  be  defined 
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that  exhibits  the  properties  of  a  Lyapunov  function  for  any  state.  If 
such  a  Lyapunov  function  can  be  defined,  the  closed  loop  system  will 
be  globally  asymptotically  stable.  To  be  a  robust  control  law,  the 
Lyapunov  function  must  also  maintain  the  properties  of  a  Lyapunov 
function  over  the  range  of  expected  system  unknowns. 

Sliding  mode  controllers  are  robust  to  both  parametric 
uncertainties  (e.g.,  imprecision  on  the  mass  properties)  and  the 

presence  of  unmodelled  dynamics  (e.g.,  structural  resonant  modes). 
System  models  are  not  required  to  be  extensive  as  long  as  the 
bounds  on  the  uncertainties  are  known.  Also,  the  less  dominant 
dynamics  are  not  required  to  be  modelled  with  high  precision. 

Instead,  they  can  be  treated  as  disturbances  as  long  as  the  upper 

bounds  on  their  magnitudes  are  known.  For  an  aerocapture  vehicle, 
for  example,  the  accelerations  due  to  aerodynamic  torques  are 
assumed  to  be  significantly  smaller  in  magnitude  than  the 

accelerations  due  to  the  control  actuators.  Being  less  dominant,  these 
aerodynamic  accelerations  can  be  treated  as  disturbances  and  an 
extensive  aerodynamic  model  is  not  required.  As  long  as  there  are 
known  bounds  on  the  system  unknowns,  a  globally  stable  sliding 
mode  control  law  can  be  designed  which  is  robust  over  the  bounded 
range  of  system  unknowns. 

3.3.1  Single-Input,  Single-Output  System  Model 

To  best  understand  the  sliding  mode  concepts  and  for  ease  of 
explanation,  first  consider  a  general,  single-input  system  [18]  with 
dynamics  described  by  the  equation 
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x(")  =  f(x)  +  b(x)  u  +  d(t) 


(3.4) 


where  the  scalar  x  is  the  output  to  be  controlled  (e.g.,  position  of  a 
body  in  one  dimension),  x(")  denotes  the  n-th  derivative  of  x,  u  is  the 
control  input  (e.g.,  imposed  force  or  acceleration),  and 
X  =  [x,  X,...,  x("-i)]'^  is  the  state  vector.  With  parametric  unknowns 
having  bounded  limits,  the  typically  nonlinear  function,  f(x),  and  the 
control  gain,  b(x),  are  not  exactly  known,  but  upper  bounds  on  the 
magnitude  of  the  imprecision  in  f(x)  and  b(x)  are  known. 
Unmodelled  dynamics  are  accounted  for  in  the  disturbance  term, 
d(t),  whose  magnitude  is  also  unknown  but  upper  bounded  by  a 
known  continuous  function  of  x  and  t.  The  control  problem,  then,  is 
to  track  a  desired  time-varying  state  ^  =  [x^,  Xd,...,  Xd^"'^^]^  despite  the 
uncertainties  in  f(x),  b(x),  and  d(t). 

As  an  example,  again  consider  the  simple  pendulum  of  the 
previous  section.  Suppose  a  perpendicular  controlling  force  u,  as 
shown  in  Figure  3.2,  is  used  to  assist  in  bringing  the  pendulum  to 
rest  at  0  =  0  (Section  3.4  presents  the  case  where  the  desired  state  is 
not  a  single  stable  equilibrium  state,  but  rather  a  time-varying 
state).  The  dynamics  of  the  system  are  now  defined  by 

mr^G -I- k9 -f- mgr  (sin  0)  -  ni  =  0  (3.5) 

or  in  the  form  of  Equation  3.4 

0  =  f(x) bu  +  d(t)  (3.6) 
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where  x  =  0  is  the  output  to  be  controlled,  the  state  vector  is 
X  =  [e  e]’',  d(t)  is  a  disturbance  term  consisting  of  any  unmodelled 
dynamics,  b  =  r,  and 


f(x)  =  -  6  -  ^  sin  0 

mr^  ^ 


(3.7) 


Figure  3.2.  Controlled  Pendulum 

Sliding  mode  control  requires  known  bounds  on  the 
imprecision  on  f(x),  b,  and  d(t).  For  this  example,  assume  the 
following  parameters  are  not  exactly  known  but  fall  within  the  given 
ranges: 
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1.0ft<r<  1.3  ft 


10!!HSi^<k<2.0  5iHSi!! 

sec  sec 

1 .0  slug  <  m  <  2.0  slug 

-0.5  -cad.  <  d  <  0.5  -taii- 

sec  2  sec  ^ 

Gravitational  acceleration  is  a  known  constant  32  ft/sec^.  Given 
these  ranges,  the  bounds  on  f(x)  can  be  calculated  from  Equation  3.7 
as 


(-  2.0  0  -  32.0  sin  g)  ^  <  f  (x)  <  (-  0.3  0  -  24.6  sin  g)  -Cad_ 

sec  2  sec  ^ 


3.3.2  Sliding  Surface 

Sliding  mode  control  design  begins  by  reducing  the  tracking 
problem  lo  a  first  order  stabilization  problem.  A  time  varying 
"sliding  surface"  is  defined  in  the  state-space  by  the  scalar  equation 

s(x,t)  =  *x  =  0  (3.8) 

where  the  scalar  x  represents  the  tracking  error  in  x  (i.e.,  x  =  x  -  x^), 
d/dt  is  the  derivative  operator,  and  X  is  a  positive  constant  with 
units  of  1/time.  Equation  3.8  suggests  that  once  on  the  sliding 
surface,  the  tracking  error  moves  exponentially  to  zero  with  a  time 
constant  (n-l)/X.  While  on  this  surface,  the  system  is  said  to  be 
behaving  in  the  "sliding  mode". 
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For  the  second-order  (i.e.,  n=2)  pendulum  example,  the  sliding 
surface  is  defined  by  the  equation 


s(x,t)  =  -I-  Xj  0  =  0 

=  0  -I-  xe  =  0 
=  0-0d  +  X(0-0d)  =  O 


For  the  simple  problem  of  bringing  the  pendulum  to  rest  at  the 
bottom  of  its  arc,  0d  =  0d  =  0  so  that  the  sliding  surface  is  given  by 

s(x,t)  =  0  +  X0  =  0  (3.10) 

With  X  a  positive  constant,  Equation  3.10  is  the  equation  of  a  line  in 
the  phase  plane  passing  through  the  origin  as  shown  in  Figure  3.3. 
The  general  solution  of  this  simple  differential  equation  is  of  the 
form 


0(t)=e„e>-' 

so  that 

0(t)  =  -0oXe-^' 

These  solutions  suggest  that  when  the  state  trajectory  is  on  the 
sliding  surface,  it  will  decay  exponentially  towards  the  desired  state 
Xd  =  [0  OF  with  a  time  constant  1/X. 

Since  the  behavior  on  the  sliding  surface  is  for  a  state 
trajectory  to  move  exponentially  towards  its  desired  state,  the  task 
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of  sliding  mode  control  is  to  drive  the  state  trajectory  to  this  sliding 
surface,  and  once  on  it,  to  maintain  it  on  the  sliding  surface.  To 
accomplish  its  task,  sliding  mode  control  attempts  to  design  a  control 
law  so  that  a  Lyapunov  function  based  on  the  closed  loop  system 
dynamics  can  be  defined  to  show  that  the  sliding  surface  is 
asymptotically  stable. 


X  i 

1 

X 

-  X,  slope 

Figure  3.3.  Sliding  Surface  Through  Origin 


As  presented  in  references  [18,19],  the  Lyapunov  function  used 
in  sliding  mode  control  is  the  positive  definite  function  s^  with  s  = 
s(x,t)  as  defined  in  Equation  3.8.  To  be  a  Lyapunov  function,  the 
derivative  of  s^  must  be  negative  despite  system  uncertainties.  This 
is  guaranteed  when  the  "sliding  condition" 

(3.11) 


is  satisfied  where  T|  is  a  positive  constant.  This  particular  inequality 
also  guarantees  that  if  initially  the  system  is  not  on  the  sliding 
surface  (i.e.,  s^^),  the  surface  will  be  reached  in  a  finite  time  less 
than 


freach  ^ 


I  s(x(0),  o)  I 

Tj 


(3.12) 


With  the  sliding  condition  (Equation  3.11)  met,  s^  is  a  Lyapunov 
function  of  the  closed  loop  dynamics  and  the  sliding  surface  is 
attractive.  Equation  3.11  ensures  all  trajectories  point  toward  the 
sliding  surface,  as  illustrated  in  Figure  3.4,  and  that  once  on  the 
surface,  the  state  trajectory  stays  on  the  surface. 
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Figure  3.4.  Attractive  Sliding  Surface 


A  typical  state  trajectory  for  a  second-order  system,  such  as 
the  pendulum,  satisfying  the  sliding  condition  (Equation  3.11)  is 
depicted  in  Figure  3.5.  From  an  initial  condition  off  the  sliding 
surface,  the  trajectory  reaches  the  surface  in  a  finite  time  less  than 
that  specified  by  Equation  3.12,  and  then  slides  exponentially  along 
the  surface  to  Xd  with  a  time  constant  of  l/A.. 


3.3.3  Switching  Sliding  Mode  Control 

The  desired  sliding  mode  control  law  for  the  general  single¬ 
input  case  is  a  switching  (i.e.,  discontinuous)  design  of  the  form 

u  =  i[u-ksgn(s)]  (3  13) 

b 

where  the  sign  function  sgn(s)  is  defined  by 

sgn(s)  =  -1  for  s  <  0 
sgn(s)  =  +1  for  s  ^  0 
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In  Equation  3.13,  u,  to  be  derived  below  for  the  pendulum,  is  the  best 
approximation  of  a  continuous  control  law  that  would  maintain  the 
system  trajectory  on  the  sliding  surface,  while  the  discontinuous 
term  ksgn(s)  satisfies  the  sliding  condition  (i.e..  Equation  3.11) 
despite  the  system  uncertainties,  thus  ensuring  state  trajectories  off 
the  sliding  surface  are  driven  back  to  the  surface.  The  remainder  of 
this  section  will  be  used  to  show  how  sliding  mode  control  design  is 
accomplished  by  designing  a  control  law  of  the  form  of  Equation  3.13 
for  the  simple  pendulum  example. 

Since  the  control  gain  b(x)  is  multiplicative  in  the  dynamics, 
the  estimated  control  gain,  b,  is  typically  chosen  as  the  geometric 
mean  of  the  bounds  on  b(x) 


b  =  Vbmin  bi 


min  «-’max 


(3.14) 


For  the  pendulum. 


b  =  V(l.0Hl.'^)  =  1.14  ft 


The  estimated  control  u  is  chosen  to  maintain  s(x,t)  =  0  which 
means  the  derivative  of  s(x,t)  must  also  equal  zero.  Differentiating 
Equation  3.10  yields 


s(x,t)  =  0  +  =  0 


Substituting  the  expression  for  the  estimated  closed  loop  dynamics  of 
the  pendulum  from  Equation  3.6  gives 


s(x,t)  =  (f  +  b  u  H  d)  +  A,0  =  0 


(3.15) 


Solving  Equation  3.15,  u  for  the  pendulum  is  found  to  be 


u  =  -  f  -  d  -  XO  (3.16) 

where  b  has  been  accounted  for  in  the  overall  control  law.  Equation 
3.13.  The  estimates  f  and  d  aie  selected  as  the  average  of  their  given 
bounds: 


f  =  ^  (-  2.0  0  -  32.0  sin  0  -  0.3  0  -  24.6  sin  0) 

2 

f  =  -  1.15  0  -28.3  sin0 

and 

d  =  l  (-  0.5  +  0.5)  =  0 

The  choice  of  X  depends  on  the  type  and  characteristics  of  the 
system  to  be  controlled.  From  reference  [19j,  for  a  mechanical 
system  X  must  be  smaller  than  all  of  the  following:  the  frequency  of 
the  lowest  unmodclled  resonant  mode,  the  sampling  frequency,  and 
the  inverse  of  the  largest  unmodelled  time  delay.  For  this  pendulum 
example,  assume  the  sampling  frequency  drives  the  strictest  criteria. 
A  general  rule  of  thumb  for  choosing  X  based  on  the  sampling 
frequency  is 


_  sampling  rate 
5 
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For  a  sampling  rate  of  100  hz,  X  =  20. 

The  gain  k  in  Equation  3.13  must  be  large  enough  to  satisfy  the 
sliding  condition.  Equation  3.11.  As  derived  in  reference  [19],  a  gain 
that  guarantees  the  sliding  condition  will  be  satisfied  can  be 
calculated  from 


k{x)  =  p(F  +  D  +  'n)  +  (|3  -  l)|u|  (3.17) 


where 


F  =  max  I  f  -  f  I 
D  =  max  I  d  -  d 


(3.18) 


are  the  maximum  deviations  of  f(x)  and  d(t)  from  their  estimates,  f 
and  d.  In  Equation  3.17,  P  is  considered  the  gain  margin  and  is 
calculated  by 


P 


(3.19) 


For  the  pendulum  example. 


P  =  fL3  =  1.14 
F  =  0.85  +  3.7  |sin  e| 
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D  =  0.5 


The  only  requirement  on  the  selection  of  tj  is  that  it  be  a  positive 
constant.  It  may  be  selected  to  provide  a  desired  reach  time 
according  to  Equation  3,12.  For  this  pendulum  example,  ti=20. 

With  the  switching  control  law  complete,  a  simulation  of  the 
pendulum  is  now  run.  To  demonstrate  the  control  law's  robust 
performance  in  the  presence  of  unknowns  and  disturbances,  the 
actual  values  of  r,  k,  m,  and  d  are  varied  over  their  full  possible 
ranges.  Releasing  the  pendulum  from  an  angle  of  45  deg,  the  control 
law.  Equation  3.13,  attempts  to  assist  the  pendulum  in  coming  to  rest 
at  0=0.  Figure  3.6  plots  the  trajectory  of  the  pendulum  vs.  time.  The 
pendulum  achieves  0=0  in  less  than  0.5  sec  and  stays  close  to  0=0. 
Figure  3.7  plots  the  value  of  the  sliding  variable  s(x,t)  vs.  time.  The 
state  trajectory  of  the  pendulum  is  initially  off  the  sliding  surface, 
s(x,0)  ^  0,  but  the  sliding  surface  is  reached  in  0.4  sec,  agreeing  with 
the  criteria  set  by  Equation  3.12: 

<  1.^^  =  IMiiM  =  2fl/45._ZL.\  =  0.8  sec 

11  T1  20l  180/ 
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Figure  3.6.  Pendulum 
(Switching  C< 


5: 


0.0 
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1.0 
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time  (sec) 

Figure  3.7.  Sliding  Variable  vs.  Time 
(Switching  Control  Law) 
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3.3.4  Continuous  Approximation  of  the  Switching  Sliding 
Mode  Control  Law 


Taking  a  closer  look  at  Figure  3.7,  Figure  3.8  shows  the 
behavior  of  s(x,t)  near  the  sliding  surface.  Since  the  second  term  of 
the  sliding  mode  control  law.  Equation  3.13,  is  discontinuous  across 
the  sliding  surface,  control  chattering  is  encountered  during  its 
implementation. 


Figure  3.8.  Sliding  Variable  vs.  Time 
(Expanded  View) 
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Such  chattering  is  undesirable  since  it  results  in  a  large  amount 
of  control  activity,  as  illustrated  in  Figure  3.9,  and  may  excite  high 
frequency  dynamics  neglected  in  the  course  of  modelling. 


time  (sec) 

Figure  3.9.  Control  Activity 
(Switching  Control  Law) 


Chattering  can  be  eliminated,  however,  with  a  continuous 
approximation  of  Equation  3,13.  A  thin  boundary  layer  about  the 
sliding  surface  is  introduced  such  that  outside  the  boundary  layer, 
the  control  law  is  as  in  Equation  3.13,  guaranteeing  the  boundary 
layer  is  attractive.  Inside  the  boundary  layer,  a  continuous  function 
is  required  in  place  of  sgn(s)  to  make  a  smooth  transition  across  the 
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sliding  surface.  These  characteristics  are  achieved  by  replacing  the 
sgn(s)  term  in  Equation  3.13  with  the  saturation  function  sat(s/d>) 
where  0  is  the  boundary  layer  thickness  and  sat(s/d>)  is  defined  by 


sat(s/<I>)  =  s/d>  for  1  s/O  1  <  1 

sat(s/<t>)  =  sgn(s/^>)  otherwise 

Figure  3.10  graphically  illustrates  the  definition  of  sat(s/<t>). 


(3.20) 


Figure  3.10.  Saturation  Function,  sat(s/<b) 


From  the  derivation  presented  in  section  7.2  of  reference  [19],  the 
sliding  condition  can  be  met  by  choosing  0(t)  as  the  solution  to  the 
appropriate  differential  equation  from 


<i)(t)  = -XO  +  p  k(ii)  ifk(xd)>^  (3  21) 
^t)  =  otherwise 

p2  p 

with  initial  condition 


<i>(o)  =  p 


k(^) 

X 


The  new,  continuous  control  law  is  then 


u  =  i-  [u  -  k‘  (x)  sat(s/0)] 
b 


(3.22) 


where  the  scalar  gain  on  the  saturation  term  is  found  from 

k'(x)  =k(x)- k(xd)  +  ^  (3  23) 

and  k(x)  and  k(Xd)  are  calculated  from  Equation  3.17.  The  continuous 
control  law  (Equation  3.22)  is  essentially  the  same  as  the  switching 
control  law  (Equation  3.13)  outside  the  boundary  layer.  Inside  the 
boundary  layer.  Equation  3.22  is  basically  a  proportional  control  law. 

Using  this  continuous  approximation  to  the  control  law  for  the 
pendulum  example.  Figure  3.11  shows  the  behavior  of  the  pendulum 
for  a  repeat  of  the  previous  simulation.  The  effort  to  eliminate 
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chattering  results  in  the  continuous  sliding  mode  control  law 
(Equation  3.22)  with  tracking  to  within  the  boundary  layer  rather 
than  the  "perfect"  tracking  of  the  switching  sliding  mode  control  law. 
Figure  3.12  shows  the  behavior  of  s(x,t)  without  chattering  near  the 
sliding  surface  while  Figure  3.13  shows  the  reduced  smooth  control 
activity  required. 


Figure  3.11.  Pendulum  Trajectory  with  0d=O 
(Continuous  Control  Law) 
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(rads/sec) 


Figure  3.12.  Sliding  Variable  vs.  Time 
(Continuous  Control  Law) 
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cont 


0  1  2 


time  (sec) 

Figure  3.13.  Control  Activity 
(Continuous  Control  Law) 
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The  state  trajectory  in  the  phase  plane  is  shown  in  Figure  3.14. 
The  line  through  the  origin  with  a  slope  of  -20  represents  the  sliding 
surface.  The  trajectory  reaches  the  sliding  surface  asymptotically 
and  then  slides  along  the  surface  to  the  origin. 


Figure  3.14.  State  Trajectory  in  Phase  Plane 
(Continuous  Control  Law) 
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3.4  The  Pendulum  Example  Continued 

In  the  previous  section,  sliding  mode  control  was  used  to  assist 
the  pendulum  in  achieving  a  stable,  equilibrium  condition  (i.e.,  x<i  = 
[0  0]'^).  The  effect  of  the  control  law  was,  essentially,  to  enhance  the 
damping  on  the  pendulum  so  that  it  reached  the  equilibrium  state 
sooner  than  it  would  without  the  control.  Rarely  does  the  desired 
state  of  a  system  correspond  to  a  single  equilibrium  state.  The  goal 
of  any  control  law  design  is  to  make  a  system  perform  in  a  desired 
behavior  whether  that  behavior  is  natural  for  the  system  or  not.  In 
this  section,  therefore,  sliding  mode  control  is  used  to  control  the 
pendulum  to  track  a  varying  trajectory. 

Though  the  desired  state,  x<i  =  [0d(t)  0d(t)]^,  does  not  represent 
a  single  stable  state  as  in  the  previous  section,  with  the  correct 
control  law,  the  time  varying  sliding  surface  does  represent  an 
asymptotically  stable  surface.  All  state  trajectories  lead  to  the 
surface,  and  once  on  the  surface,  the  trajectory  slides  exponentially 
along  the  sliding  surface  to  the  desired  state. 

For  this  example,  the  pendulum  is  started  at  rest  at  the  bottom 
of  its  arc,  x(0)  =  [0  0]^^.  The  desired  trajectory  is  given  by  the 

equation 


(3.24) 


so  that  the  pendulum  must  track  a  slow  swing  from  0d(O)  =  +45  deg 
to  0d(5  sec)  =  -45  deg  and  back  to  0d(lO  sec)  =  +45  deg. 


The  desired  control  law  is  still  of  the  form  of  Equation  3.22, 
with  the  estimates  b,  f,  d  and  the  parameters  (3,  F,  D,  X,  and  q  the 
same  as  before.  The  sliding  surface  is  defined  by 


s(x,t)  =  e  +  Xe  =  0  (3.25) 

as  before,  but  since  is  now  time  varying,  the  sliding  surface  no 
longer  passes  through  the  origin  of  the  phase  plane,  but  rather  is  a 
time  varying  line  through  Xd  with  slope  -X  =  -20. 

The  best  approximation  u  of  a  continuous  control  law  that 
would  maintain  the  state  trajectory  on  the  sliding  surface  is  again 
found  by  requiring  the  derivative  of  s(x,t)  to  be  zero.  Differentiating 
Equation  3.25, 


s(x,t)  =  0  -I-  X0  =  0 

=  e-9<,  +  x(9-y 


Substituting  the  equation  for  the  estimated  closed  loop  dynamics  of 
the  pendulum  for  0, 


s(x,t)  =  (f  b  u  +  d)  -  0d  +  >^0  -  0d)  (3.26) 

Solving  Equation  3.26  for  u  results  in 

u  =  -  f  -  d -I- 9d- X(0  -  0d)  (3.27) 
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with  b  accounted  for  in  the  overall  control  law.  Differentiating 
Equation  3.24  twice,  the  desired  angular  velocity  and  acceleration  are 


(3.28) 


and  all  terms  in  Equation  3.27  are  known.  With  the  second  term  of 
Equation  3.22  calculated  as  in  the  previous  section,  the  continuous 
sliding  mode  control  law  design  is  complete.  Figure  3.15  shows  the 
tracking  response  of  the  pendulum  using  this  sliding  mode  control 
law. 


time  (sec) 

Figure  3.15.  Pendulum  Trajectory  with  Varying  Gj 
(Continuous  Control  Law) 
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Figure  3.16  shows  the  sliding  surface  was  reached  in  0.7  sec,  again  in 
agreement  with  Equation  3.12,  with  tracking  to  within  the  boundary 
layer  after  reaching  the  sliding  surface. 


time  (sec) 

Figure  3.16.  Sliding  Variable  and 
Boundary  Layer  vs.  Time 

3.5  Application  of  Sliding  Mode  to  Attitude  Control 
for  Aerocapture 

Sliding  mode  techniques  can  be  applied  directly  to  multi-input 
systems  with  differences  from  the  single-input  case  due  only  to  the 
increased  mathematical  complexity.  The  general  form  of  a  multi¬ 
input,  multi-output  sliding  mode  control  law  is  presented  in  Section 
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3.5.1.  A  system  model  for  the  AFE  vehicle  in  a  form  suitable  for 
application  of  sliding  mode  techniques  is  developed  in  Section  3.5.2 
along  with  the  design  of  the  desired  control  law. 

3.5.1  Multi-Input,  Multi-Output  Sliding  Mode  Control 

A  general  multi-input,  multi-output  dynamical  system  can  be 
modelled  by  the  system  of  equations 

m 

=  fi  (x)  +  ^  bij  (x)  Uj  +  di (t)  (3 .2r , 

j  =  i 
i  =  1,-  • 

where  the  state  vector  x  is  composed  of  the  xi’s  and  their  first  (ni-1) 
derivatives,  di(t)  is  the  disturbance  component  associated  with  xi,  u 
is  the  control  input  vector  of  components  uj,  and  all  the  bij's  form  the 
input  matrix  B .  Uncertainty  on  the  elements  of  the  input  matrix  are 
indicated  by  writing  B  in  the  form 

B=(i  +  A)B  (3.30) 

where  1.  is  the  m  x  m  identity  matrix,  B  is  the  estimated  input 

matrix,  and  the  matrix  A  is  a  measure  of  the  uncertainty. 

Following  the  logic  presented  in  Section  3.3  for  the 

development  of  the  control  law  for  the  single-input  case,  a 

discontinuous  switching  control  law  is  first  sought  for  the  multi¬ 
input,  multi-output  system  that  verifies  the  sliding  conditions 


(3.31) 


<-Tii|sii.  'ni>o 

z  at 

where  each  si  is  a  component  of  a  vector  s  and  is  defined  by 

Si(x,t)  =  (^  +  ^i)"‘  ^  Xi  (3.32) 

similar  to  Equation  3.8  for  single-input  systems.  By  allowing  ksgn(s) 
to  be  the  vector  of  components  kisgn(si),  the  desired  switching 
control  law  can  be  written  in  the  form 

u  =  B  ^  [u  -  k  sgn(s)]  (3.33) 

where,  similar  to  u  in  Equation  3.13,  u  represents  the  best 
approximation  of  a  continuous  control  law  that  would  maintain  the 
sliding  mode  for  each  dimension.  The  gains  k  are  found  by  solving 
the  system  of  equations 

(1  -  Dij)  ki  +  X  Dij  kj  =  Fi  +  X  D«  I  x^"'>  -  fj  I  +  Tii  +  di„„  (3.34) 

J=1 

where  each  Fi  can  be  calculated  from  Equation  3.18,  dim  ax  is  the 
magnitude  of  the  maximum  disturbance,  the  elements  Dij  are 
determined  from 

Dij  =  max|Aij|  (3.35) 
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and 


ri  di 


(3.36) 


As  for  the  single-input  case,  a  continuous  approximation  to  the 
control  law  (Equation  3.33)  is  required  to  eliminate  chattering  upon 
implementation  of  the  multi-input,  multi-output  control  law.  The 
desired  form  is 


u  =  B 


-1 


u  -  k  sa 


©I 


(3.37) 


where  the  term  k  sat(s/  O)  represents  a  vector  containing  the 
components  kisat(si/Oi)  with  the  saturation  function  defined  as  in 
Section  3.3  and  each  boundary  layer  thickness  0 i  determined  by  the 
desired  tracking  precision.  The  aerocapture  attitude  control  law 
design  outlined  below  presents  one  method  of  determining  the  Oi's. 

3.5.2  Aerocapture  Attitude  Control 

The  first  step  in  designing  a  control  law  of  the  form  of  Equation 
3.37  for  a  multi-input,  multi-output  system  is  to  develop  a  dynamics 
model  in  the  form  of  Equation  3.29  with  the  input  matrix  B  defined 
by  Equation  3.30.  Once  an  appropriate  system  model  is  developed 
and  bounds  on  the  unknowns  are  determined,  the  equations 
presented  in  the  previous  section  can  be  used  directly  to  complete 
the  control  law  design. 
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Ignoring  aerodynamic  interaction,  the  dynamics  of  an 
aerocapture  vehicle  like  the  one  to  be  used  in  the  AFE  can  be 
described  by 


t  =  dli  =  ?li+ wxH 

dt  3t 


(3.38) 


where  T,  H,  and  (0  represent  the  three-dimensional  torque  due  to  the 
actuators,  the  angular  momentum,  and  the  angular  velocity  vector, 
respectively,  with  the  angular  velocities  measured  with  respect  to 
the  body  referenced  system.  Appendix  A  compares  the  body 
referenced  system  to  the  velocity  vector  referenced  system  (i.e., 
velocity  angles).  Angular  momentum  about  the  center  of  mass  can 
be  expressed  by  H=IoJ  [30]  where  I,  the  inertia  matrix  of  the  vehicle, 
is  of  the  form 

Ixx  Ixy 
I  —  Ixy  lyy 
Ixz  lyz 

Assuming  I  changes  slowly  enough  to  be  considered  constant  over 
the  interval  of  interest. 


^xz 


lyz 


(3.39) 


3(lm)  _  I  aw 
at 


and  Equation  3.38  can  be  rewritten  as 


The  goal  of  the  attitude  control  system  is  to  control  the 
vehicle's  velocity  angles  which  are  used  to  define  the  vector 

e=[<l>  a  p]'^ 


so  that  CO  =  R  0  and  cb  =  R  0 

where 

cos  P  cos  a 

-sin  P 

cos  P  sin  a 

R  = 

sin  P  cos  a 

cos  P 

sin  P  sin  a 

-sin  a 

0 

cos  a 

is  the  rotation  matrix,  calculated  in  Appendix  A,  that  is  required  to 
convert  from  a  body  referenced  system  to  a  velocity  angle 
referenced  system.  Substituting  the  expressions  for  (0  and  cb  into 
Equation  3.40  and  rearranging  yields  the  second-order  system 

0  =  -R  e  X  iR**  e)-!-  R  r’  Tin  +  d(t)  0.42) 

where  d(t)  has  been  added  to  represent  the  disturbances  due  to 
unmodelled  dynamics  and  Tin  is  the  torque  due  to  the  control  inputs 
(i.e.,  jet  firings).  With  the  acceleration  due  to  dynamic  coupling 
already  included  in  the  first  term  of  Equation  3.42,  the  acceleration 
due  to  control  inputs  can  be  assumed  to  be  u  =  I  ^  Tin,  which 
represents  the  desired  control  (i.e.,  the  body  referenced  angular 
acceleration  commands  to  be  used  by  the  actuator  selection 


algorithm).  Rearranging,  Equation  3.42  can  be  written  in  the  desired 
form 


0  =  f  +  Bu  +  d 


(3.43) 


where 


f  =  -RrV^exiR'^e)  (3.44) 

and 

Bu  =  R  r*  Tin 

From  Equation  3.30,  the  desired  form  of  B  is  B  =  (1  +  A)B  with  1  the 
3x3  identity  matrix.  The  rotation  matrix  R  is  strictly  a  function  of 
the  velocity  angles.  Assuming  accurate  knowledge  of  these  angles,  R 
can  be  assumed  known.  The  mass  properties,  however,  are 
unknown,  and  the  actual  inertia  matrix  is  written  as  I  =  I  +  51  where 
I  is  the  current  estimate  and  51  represents  the  error  or  difference 
between  the  actual  and  estimated  inertia  matrix.  B  can  then  be  put 
in  the  desired  form  with 


5  =  1  (3.45) 

a  =  r(i  +  si)'‘ I -  i 


With  an  appropriate  system  model  established  as  Equation 
3.43,  a  discontinuous  3-dimensional  control  law  for  this  second-order 
system  of  the  form  of  Equation  3.33  can  be  designed.  Following  the 


single-input,  second-order  example  presented  in  Section  3.4,  u  is 
found  to  be 


^ 

U  =  0r - f  -  d 


(3.46) 


with 


er  =  0d'Ae  (3.47) 

where  A  is  a  diagonal  matrix  containing  the  Xj's  from  Equation  3.32, 
and  0d  represents  the  desired  velocity  angle  accelerations.  All  that 
remains  to  finish  the  switching  control  iaw  design  is  to  calculate 
estimates  and  bounds  for  I,B,  and  d(t)  based  on  the  expected 
atmospheric  densities  and  vehicle  properties  for  the  given 
aerocapture  mission  and  vehicle.  Calculations  for  these  estimates 
and  bounds  for  an  Earth  aerocapture  using  the  AFE  vehicle  are 
presented  in  Appendix  B. 

The  final  attitude  control  law  design  in  the  continuous  form  of 
Equation  3.37  can  now  be  completed  by  selecting  the  thickness  of  the 
boundary  layers.  The  elements  of  the  vector 

represent  the  boundary  layer  thickness  about  the  sliding  surfaces 
corresponding  to  each  of  the  velocity  angles.  For  simplicity,  a 
constant  thickness  is  chosen  for  each  Oi.  For  this  design,  the  choice 
of  each  <I>i  is  based  on  the  desire  to  track  the  velocity  angles  to 
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within  2  degs  of  the  commanded  angles.  Since  this  is  a  second-order 
system,  from  Equation  3.32 


Si(x,t)  =  0i  +  XiGi 


As  in  the  example  in  Section  3.3,  the  sampling  rate  is  assumed  to 
provide  the  strictest  criteria  on  the  choice  of  each  X.,.  With  a 
sampling  rate  of  25  hz,  Xi=  5.  Assuming  0i  is  negligible  with  respect 
to  Xi0i  and  angular  tracking  is  desired  to  be  within  2  deg 
(|0i|<O.O35  rads),  after  each  sliding  surface  is  reached  si  should  be 
within  the  range  I  Si  1^0.17,  The  thickness  of  each  boundary,  Oj,  is 
therefore  selected  to  be  0.17,  or 


0.17 

0.17 

-0.17 


rads. 

sec 


With  selection  of  O,  design  of  the  first  element  of  the  aerocapture 
attitude  control  system,  the  control  law,  is  complete. 

The  same  steps  described  above  can  be  used  to  design  the 
control  law  for  any  aerocapture  mission  or  vehicle.  The  only 
differences  from  one  mission  or  vehicle  to  the  next  are  the  estimates 
and  bounds  on  I,  B,  and  d(t)  and  the  choice  of  each  Xi  and  Oj. 
Estimates  and  bounds  on  I,B,  and  d(t)  are  dependent  on  the 
atmosphere  and  vehicle  properties  for  a  given  mission,  while  Xi  and 
<Di  are  based  on  the  system  cycling  rates  and  desired  tracking 
response. 
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Chapter  4 

Actuator  Selection 


4.1  Introduction 

The  task  of  the  actuator  selection  algorithm  is  to  determine 
which  actuators  to  employ,  and  to  what  extent,  in  order  to  generate 
the  angular  accelerations  commanded  by  the  control  law.  In  order 
for  the  algorithm  to  make  the  proper  selection,  accurate  knowledge 
of  the  effectiveness  of  each  actuator  is  required.  The  selection 
algorithm  must  also  be  able  to  adapt  to  updated  knowledge  of  those 
parameters  which  affect  the  control  authority  of  the  actuators  (e.g., 
vehicle  mass  properties). 

Section  4.2  presents  the  general  actuator  selection  problem  as  a 
linear  system  of  equations  whose  solution  is  the  optimal  mix  of 
actuator  activity  that  generates  the  commanded  angular 
accelerations.  Section  4.3  describes  how  linear  programming  is  used 
with  the  simplex  method  to  solve  the  actuator  selection  problem. 

4.2  Actuator  Selection  Problem  as  a  Linear  System 
of  Equations 

The  actuator  selection  problem  can  be  represented  as  the 
problem  of  finding  a  solution  to  a  system  of  linear  equations  equal  to 
the  requested  angular  accelerations.  In  addition,  it  may  be  desirable 
to  find  the  solution  that  minimizes  a  linear  cost  function 
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corresponding  to  fuel  usage.  Such  a  problem  can  be  summarized  as 
the  minimization  of  the  cost  function 

n 

Z  =  ^  Cj  Xj  (4.1) 

j  =  l 

subject  to 

n 

^  ajXj  =  AR  (4.2a) 

j  =  i 

0^Xj<u/  (4,2b) 

where 

n  =  #  of  actuators  available 
cj  =  cost  factor  associated  with  use  of  actuator  #j 
xj  =  decision  variable  indicating  amount  of  action  selected 
for  actuator  #j 

.1+ 

“J  =  upper  bound  on  the  decision  variable  for  actuator  #j 
aj  =  activity  vector  specifying  control  authority  of 
actuator  #j 

AR  =  [Pc  qc  Tc  I'T  =  commanded  angular  accelerations. 

This  formulation  of  the  actuator  selection  problem  can  be 
applied  to  more  than  one  type  of  actuator.  Each  aj  represents  the 
acceleration  provided  by  actuator  #j  per  unit  amount  of  its 
employment.  The  solution  x  to  the  linear  program  is  the  vf^ctor 
containing  the  decision  variables  xj  specifying  the  amount  of 
corresponding  actuator  action  required  to  generate  the  commanded 
acceleration.  For  aerosurface  control  actuators,  decision  variables 
may  represent  aerosurface  deflection  angles.  For  RCS  jets,  the 


80 


decision  variables  represent  the  duty  cycle  (i.e.,  the  ratio  of  jet  on 
time  to  off  time),  ranging  from  0  to  1,  defining  the  fraction  of 
maximum  acceleration  required  of  a  particular  jet.  The  cost  cj  of 
using  actuator  #j  is  defined  by  the  user  and  can  vary  greatly 
depending  on  the  type  of  actuator  and  the  mission  profile.  The  cost 
of  employing  an  aerosurface  relative  to  the  cost  of  an  RCS  jet  may  be 
higher  at  high  altitudes  where  aerosurfaces  are  less  effective  than 
jets,  but  may  be  lower  at  the  lower  attitudes  where  they  are 
effective  and  reduce  the  fuel  expenditure  due  to  firing  jets. 

For  the  AFE-like  vehicle  used  in  this  thesis,  the  only  actuators 
used  for  attitude  control  are  RCS  jets.  When  new  angular 
accelerations  are  commanded  by  the  control  law,  a  new  optimization 
problem  of  the  form  of  Equations  4.1  and  4.2  must  be  solved  to 
determine  the  optimal  jet  duty  cycles.  The  scheduling  of  jet  firings  is 
then  determined  based  on  the  chosen  duty  cycles  as  will  be 
discussed  in  Chapter  6.  To  perform  this  jet  selection  task,  an  activity 
vector  aj  for  each  jet  must  first  be  calculated. 

The  activity  vector  of  each  jet  is  determined  from  the  equation 

aj  =  rMrjxfj)  (4.3) 

where  I  is  the  vehicle  inertia  matrix,  rj  is  the  position  of  jet  #j  witn 
respect  to  the  center  of  mass  of  the  vehicle,  and  Tj  is  the  thrust  of  jet 
#j.  Figure  4.1  is  a  rear  view  of  the  AFB  vehicle  used  in  this  thesis 
diagraming  the  positioning  of  the  RCS  jets.  Jets  #1  -  #8  all  provide 
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Jet  # 

Direction  of 
thrust 

1 

+Y 

2 

-  Y 

3 

+Z 

4 

-Z 

5  -  16 

+X 

Figure  4.1.  Aft  View  of  AFE  Vehicle 
(+X  axis  into  page) 


30  Ibf  of  thrust  while  jets  #9  -  #16  provide  125  Ibf  of  thrust  [31]. 
Assuming  known  position  and  thrust  of  each  jet,  accurate  estimates 
of  the  vehicle  inertia  and  center  of  mass  location  are  required  to 
properly  calculate  jet  activity  vectors.  Estimation  of  these  mass 
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properties  is  accomplished  outside  of  the  actuator  selection  algorithm 
(i.e.,  by  the  mass  property  identification  algorithm  discussed  in 
Chapter  5).  Whenever  updated  estimates  are  provided,  the  activity 
vectors  are  recalculated. 

All  that  remains  to  complete  the  formulation  of  the  jet  selection 
problem  is  the  definition  of  the  costs  associated  with  the  use  of  each 
jet.  For  the  candidate  control  system,  each  cj  is  assigned  as  a 
constant  representing  the  rate  of  fuel  consumption  during  firing  of 
jet  #j.  Assuming  a  constant  ratio  between  the  amount  of  fuel  used 
and  the  amount  of  thrust  provided,  the  125  Ibf  jets  are  assigned  a 
cost  of  4.2  while  the  30  Ibf  jets  are  assigned  a  cost  of  1. 

Several  algorithms  surveyed  in  reference  [28]  can  be  used  to 
solve  the  system  represented  by  Equation  4.2.  Of  these,  only  the 
linear  programming  approach  solves  Equation  4.2  while  minimizing 
Equation  4.1.  Because  of  this,  linear  programming  was  chosen  for  the 
jet  selection  element  of  the  candidate  aerocapture  attitude  control 
system. 

4.3  Linear  Programming  and  the  Simplex  Method 

The  use  of  linear  progra  nr.jing  has  been  shown  [32]  to  perform 
the  task  of  adaptable  jet  selection  quite  well.  It  is  especially  well 
suited  for  hybrid  actuator  selection  for  which  the  optimal  blend  of 
different  types  of  actuators  (e.g.,  RCS  jets  and  aerosurfaces)  is  chosen 
from  a  common  pool  of  actuators  [20,21].  RCS  jets  are  the  only 
attitude  control  actuators  available  to  the  actual  AFE  vehicle.  For  the 
current  effort,  they  also  represent  the  only  type  of  actuators 
available  to  the  candidate  design.  As  an  actuator  selection  algorithm 
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for  jets  only,  linear  programming  is  still  attractive,  having  been 
shown  in  flight  test  to  select  the  most  fuel  efficient  mix  of  jet  firings 
that  implements  the  commanded  angular  accelerations.  Also,  by 
employing  linear  programming  for  jet  selection  for  this  example,  this 
attitude  control  system  can  easily  be  used  for  control  of  a  vehicle 
with  varied  actuators  without  changing  the  system  structure. 

The  combination  of  Equations  4.1  and  4.2  is  the  formulation  of 
a  general  linear  programming  problem.  To  solve  this  problem,  an 
upper  bounding  simplex  method  is  employed.  Section  4.3.1  outlines 
the  basic  simplex  method  for  a  general  linear  programming  problem. 
Section  4.3.2  describes  how  the  simplex  method  is  altered  v/hen  the 
decision  variables  (i.e.,  amounts  of  actuator  activity)  are  upper 
bounded  as  in  the  actuator  selection  problem. 

4.3.1  Basic  Linear  Programming  Problem  and  Simplex 
Method 

Like  the  actuator  selection  problem  presented  in  Section  4.2, 
the  typical  linear  programming  problem  is  to  find  the  one  solution  to 
an  underdetermined  system  of  linear  equations  having  many 
possible  solutions  which  minimizes  a  linear  cost  function.  In 
equation  form  [33],  the  problem  is  to  minimize  z  -  c  ^  x  subject  to 

Ax=b,Xj>0  (4.4) 


where  x  is  the  solution  vector  having  n  components,  c  is  a  vector 
having  n  cost  components,  z  =  cixi  +  ...  +  c^x^  is  the  cost  function  to 
be  minimized,  A  is  an  m  by  n  matrix  with  m  ^n,  and  b  is  a  column 
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vector  of  m  components.  Comparing  this  formulation  to  Equations 
4.1  and  4.2,  the  columns  of  the  matrix  A  are  the  activity  vectors 
associated  with  each  decision  variable  xj,  and  b  is  the  requested 
output  AR .  The  main  difference  between  the  two  formulations  is 
that  the  xj's  are  not  upper  bounded  in  Equation  4,4  as  they  are  in 
Equation  4.2. 

An  important  property  of  the  solution  vector  x  that  minimizes 
the  cost  function  is  that  at  least  n  -  m  of  its  components  are  zero 
[33].  To  illustrate  this  property,  consider  the  simple  problem  of 
minimizing  the  cost  function  z  =  2xi  +  X2  subject  to  the  equation  xi  + 
2x2  =  4.  Written  in  the  form  of  Equation  4.4,  the  problem  is: 


minimize 


subject  to 


z  =  2x1  +  X2 


[12] 


xi 

X2 


=  4 


(4.5a) 


(4.5b) 


with  m  =  1,  «  =  2  (m  <  n).  In  Figure  4.2,  the  line  segment  PQ 
represents  the  "feasible  set"  of  all  possible  solutions  to  the  equality 
constraint  (Equation  4.5b)  with  the  additional  constraint  >f  xj  >  0.  For 
a  cost  of  z=0.  Figure  4,2  shows  that  the  line  of  constant  cost  does  not 
intersect  the  feasible  set.  The  first  intersection  occurs  at  point  P 
when  the  cost  is  increased  to  2.  This  intersection  point, 
corresponding  to  x  =[0  2]T,  represents  the  solution  to  the  equality 
constraints  that  minimizes  the  cost  function  and,  therefore,  the 
desired  solution  to  the  linear  program.  The  solution  also  agrees  with 
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the  property  of  having  n  -  m  =  1  components  equal  to  zero.  This 
property  holds  for  problems  of  higher  dimension. 


Figure  4.2.  Feasible  Set  and  Cost  Curves 


As  outlined  in  reference  [33],  the  simplex  method  can  be  used 
to  solve  the  linear  programming  problem.  Given  the  property 
discussed  above,  the  task  of  the  simplex  method  is  to  determine 
where  the  zeros  belong  in  the  solution.  The  feasible  set  of  solutions 
can  be  imagined  to  be  a  polyhedron  of  many  edges  and  faces  with 
the  corners  representing  a  solution  to  the  equality  constraints  in 
Equation  4.4  with  only  m  nonzero  components.  The  basic  approach 
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of  simplex  is  to  start  at  an  extreme  point  or  comer  of  the  feasible  set. 
Since  this  starting  point  probably  does  not  minimize  the  cost 
function,  the  next  step  in  simplex  is  to  move  to  a  neighboring  corner 
of  the  feasible  set  where  the  cost  is  lower.  To  find  a  neighboring 
corner,  one  of  the  zero  components  of  x  is  allowed  to  become  positive 
while  the  remaining  zero  components  are  unchanged.  As  the 
selected  component  is  increased,  the  original  nonzero  components 
must  adjust  in  order  for  the  equality  constraints  in  Equation  4.4  to 
still  hold.  During  the  adjustment,  some  of  these  components  may 
decrease  toward  zero.  Once  the  new  component  of  x  is  increased  to 
the  point  where  one  of  the  original  nonzero  components  becomes 
zero,  a  neighboring  corner  has  been  found  and  the  cost  can  be 
evaluated  to  determine  if  it  has  been  reduced.  Since  there  is  a  finite 
number  of  corners  to  a  feasible  set,  a  smaller  finite  number  of  steps 
is  required  until  the  cost  is  reduced  to  its  minimum  and  the  optimal 
solution  is  found. 

As  an  example  illustrating  the  basic  approach  of  simplex, 
consider  the  following  simple  minimization  problem: 


minimize 
subject  to 


z  =  [3  4  1]  X 

[  1  1  1  ]  X  =1  ,  Xi  >  0 


(4.6a) 


(4.6b) 


The  shaded  region  in  Figure  4.3  represents  the  feasible  sei  of 
solutions  to  the  equality  constraints  in  Equation  4.6b. 
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Figure  4.3.  3-D  Feasible  Set  for  Example  Problem 

Point  P  represents  a  corner  of  the  feasible  set  with  m  =  1  positive 
components  and  n  -  m  =  2  zero  components.  Beginning  at  P  (i.e.,  x  = 
[1  0  0]T),  simplex  must  find  a  neighboring  corner  with  a  cost  lower 
than  z(P)  =  3.  The  decision  as  to  which  zero  component  to  allow  to 
become  positive  must  be  made.  If  X2  is  increased  while  keeping 
X3=0,  then  Equation  4.6b  requires 

XI  =  1  -  X2  (4.7a) 

and  the  cost  function  is 

z=3(l  -  X2)  +  4x2  =  3  +  X2  (4.7b) 
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so  that  the  cost  is  increased  as  X2  is  increased.  If,  however,  X3  is 
increased  while  X2  is  held  at  zero.  Equation  4,6b  requires 

XI  =  1  -  X3  (4.8a) 

and  the  cost  function  is 

z=3(l  -  X3)  +  X3  =  3  -  2x3  (4.8b) 

resulting  in  the  cost  being  reduced  with  the  increase  of  X3.  The 
decision  to  increase  X3  instead  of  X2  is  obvious.  Equation  4.8a  is  then 
used  to  determine  how  much  to  increase  X3  in  order  to  arrive  at  a 
comer  of  the  feasible  set.  If  X3  is  increased  to  1,  xi  is  reduced  to  zero 
and  the  point  x  =  [0  0  1]^,  labelled  R  in  Figure  4.3,  is  the 
neighboring  corner  selected  by  simplex.  Since  the  only  corners 
neighboring  R  represent  higher  costs,  point  R  is  the  optimal  solution 
to  the  linear  programming  problem  presented  by  Equation  4.6. 

4.3.2  Upper  Bounding  Simplex  for  Actuator  Selection 

Since  there  are  upper  limits  to  the  amount  a  control  actuator 
can  be  employed  (i.e,,  an  aerosurface  can  only  be  deflected  a  limited 
number  of  degrees,  and  100%  is  the  maximum  duty  cycle  of  an  RCS 
jet)  the  decisions  variables  Xj  are  upper  bounded,  and  simplex  must 
not  select  decision  variables  above  their  limits.  A  modified,  upper 
bounding  simplex  method  is,  therefore,  required  [21].  As  with  the 
basic  simplex  method,  this  modified  method  begins  with  a 
nonoptimal  solution  to  the  linear  program  with  m  nonzero 
components  of  x  representing  the  use  of  m  actuators  for  the  solution. 
An  w  by  m  matrix  B,  whose  columns  consist  of  the  activity  vectors  aj 
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associated  with  the  m  chosen  actuators,  is  considered  the  basis.  For 
three  dimensional  attitude  control,  m  =3  and  upper  bounding  simplex 
begins  with  a  3-jet  solution  with  the  basis  B  a  3  by  3  matrix 
containing  the  activity  vectors  of  the  jets  chosen  for  the  solution. 

Simplex  then  considers  one  of  the  n  -  m  unused  actuators  and 
increases  its  associated  decision  variable  until  one  of  three 
possibilities  occurs.  The  first  possibility  is  that  one  of  the  original 
decision  variables  will  go  to  zero.  When  this  happens,  if  the  cost  is 
lowered,  the  activity  vector  associated  with  the  new  nonzero  decision 
variable  is  "pivoted"  into  the  basis  in  place  of  the  activity  vector 
associated  with  the  decision  variable  driven  to  zero.  This  occurrence 
is  the  equivalent  of  a  single  step  of  the  general  simplex  method. 

The  second  possibility  during  consideration  of  the  new  actuator 
is  that  either  one  of  the  original  nonzero  decision  variables  or  the 
new  nonzero  decision  variable  reaches  its  upper  bound  before  one  of 
the  decision  variables  is  driven  to  zero.  The  resulting  solution  for 
consideration  then  consists  of  m  decision  variables  at  less  than  their 
upper  bounds  with  their  associated  activity  vectors  making  up  the 
new  basis  B ,  and  one  (or  more  after  consideration  of  several 
actuators)  decision  variable  at  its  upper  bound.  In  equation  form 

P 

BxB+X2kUk  =  AR  (4-9) 

k=  I 

where 

B  =  m  by  m  basis  matrix  consisting  of  the  m  activity 
vectors  with  decision  variables  that  are  nonzero 
and  below  their  upper  bounds 
Xb  =  vector  of  the  m  decision  variables  whose  associated 
activity  vectors  are  in  the  basis 
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p  =  total  number  of  decision  variables  at  their  upper 
bounds 

ak  =  activity  vector  with  associated  decision  variable  Uk 
at  its  upper  bound 

AR  =  [Pc  Qc  Tc  1^  =  commanded  angular  accelerations. 

If  this  resulting  solution  provides  a  lower  overall  cost,  the  solution  is 
saved  and  the  modified  simplex  method  selects  another  actuator  for 
consideration. 

The  third  possibility  can  occur  when  the  currently  saved 
solution  already  has  decision  variables  at  their  upper  bounds  as  in 
Equation  4.9.  As  a  new  actuator  is  considered,  the  increase  of  its 
associated  decision  variable  may  drive  a  decision  variable  at  its 
upper  bound  toward  zero.  The  new  decision  variable  is  increased 
until  a  new  solution  of  the  form  of  Equation  4.9  is  found.  Decision 
variables  at  t^eir  upper  bounds  may  drop  out  of  the  solution,  they 
may  become  a  decision  variable  whose  activity  vector  is  in  the  basis, 
or  decision  variables  associated  with  basis  activity  vectors  may  be 
driven  to  their  upper  bounds.  All  actuators  are  considered  until  the 
optimal  solution  of  the  form  of  Equation  4.9  is  found. 
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Chapter  5 

Mass  Property  Identification 


5.1  Introduction 

From  Chapter  4  it  is  obvious  that  the  linear  programming 
algorithm  used  for  jet  selection  requires  accurate  mass  property 
estimates  to  properly  choose  which  jets  to  fire.  In  particular,  activity 
vectors  specifying  the  acceleration  provided  by  each  jet  depend 
directly  on  accurate  knowledge  of  the  vehicle  inertia  matrix,  the 
location  of  the  jet  with  respect  to  the  vehicle  center  of  mass,  and  the 
thrust  vector  (i.e.,  thrust  level  and  direction)  for  each  jet.  For  the 
purposes  of  this  thesis,  the  location  on  the  vehicle  and  the  thrust 
vector  for  each  jet  are  assumed  known,  while  the  vehicle  inertia 
matrix  and  location  of  the  center  of  mass  must  be  estimated. 

For  a  mission  like  the  AFE  in  which  mass  properties  remain 
nearly  constant  and  accurate  estimates  are  known  a  priori,  the  initial 
estimates  are  valid  throughout  the  mission  and  updated  estimates 
are  not  required.  For  the  more  general  case  in  which  mass 
properties  may  be  unknown  a  priori  or  may  change  significantly 
during  the  mission,  accurate  and  timely  estimates  are  essential.  For 
many  potential  missions,  updates  of  mass  property  estimates  from 
the  control  centers  on  the  ground  are  not  feasible.  For  an 
aerocapture  mission  at  Mars,  for  example,  the  time  required  to 
transmit  data  from  Earth  makes  timely  updates  of  mass  property 
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estimates  impossible.  An  in-flight  mass  property  estimation 
capability  is,  therefore,  required. 

One  approach  to  updating  mass  property  estimates  is  to  begin 
with  initial  accurate  measurements  of  the  desired  mass  properties  of 
both  the  vehicle  and  the  payload.  Using  onboard  fuel  gauges  to 
measure  fuel  expenditure,  updated  estimates  can  then  be  generated 
by  monitoring  the  mass  of  the  fuel  consumed  along  with  the  fuel's 
location  in  the  vehicle.  Since  accurate  knowledge  of  the  payload 
mass  properties  is  required,  this  approach  is  not  feasible  for  the 
example  system  of  this  thesis  in  which  unknown  payload  mass 
properties  are  possible.  An  estimation  algorithm  that  doesn't  require 
prior  mass  property  knowledge  is  needed. 

Section  5.2  describes  a  nonlinear  filter  design  that  is  used  to 
generate  estimates  of  the  vehicle  inertia  matrix  and  location  of 
center  of  mass,  without  prior  knowledge  of  these  properties,  based 
on  noisy  measurements  of  the  vehicle  angular  rate  changes  due  to  jet 
firings.  A  strategy  to  select  the  jet  firings  that  provide  the  most 
information  to  the  filter  during  times  dedicated  to  mass  property 
identification  is  outlined  in  Section  5.3. 

5.2  Estimation  of  Vehicle  Inertia  Matrix  and  Center 
of  Mass  Location 

A  standard  Kalman  filter  is  used  to  estimate  system  unknowns 
based  on  a  linearized  model  of  the  system,  while  an  extended  Kalman 
filter  is  used  for  estimation  when  the  more  general  nonlinear  system 
model  must  be  retained  [34].  The  model  of  the  dynamics  of  a 
spacecraft  (i.e..  Equation  3.33)  is  nonlinear.  As  will  be  shown  below. 
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the  full  nonlinear  model  must  be  used  in  order  to  make  estimates  of 
all  the  unknown  mass  properties  (i.e.,  inertia  matrix  and  location  of 
center  of  mass).  The  mass  property  estimator  is,  therefore,  a 
nonlinear  filter  design  resembling  a  second-order,  extended  Kalman 
filter  and  is  of  the  form  first  presented  in  references  [22,23,24].  A 
model  of  the  rigid  spacecraft  dynamics  is  used  with  the  current  mass 
property  estimates  to  predict  the  output  of  the  spacecraft's  rate 
gyros  due  to  jet  firings.  The  filter  compares  these  predictions  with 
the  measured  values  and  makes  revisions  to  its  mass  property 
estimates  based  on  this  comparison.  The  following  outlines  the  filter 
design  of  [23,24]. 

In  Chapter  3,  a  model  of  the  dynamics  of  a  rigid  aerocapture 
spacecraft  was  developed  and  expressed  in  terms  of  the  velocity 
angles  and  rates.  A  similar  model  in  terms  of  the  body  referenced 
angular  rates  can  be  expressed  by  the  equation 


I(b  =  -  (co  X  Iw)  -I-  T  -I-  d 


(5.1) 


where 

I  =  the  vehicle  inertia  matrix 
X  =  applied  torques 
d  =  the  sum  of  all  disturbance  torques 
(b  =  angular  accelerations 
(0  =  angular  rates 


Assuming  mass  properties  and  jet  thrust  are  constant  over  the 
interval  At  and  that  to  =  Aco/At  where  Aco  is  the  change  in  angular 
rate.  Equation  5.1  can  be  rewritten  as 
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A(0  =  r'  -  (co  X  Iw)  +  (fj  X  Tj)  +  d]  At  (5  2) 

j  J 

where  Tj  is  the  known  thrust  of  jet  #j,  rj  is  the  displacement  of  jet  #j 
fiom  the  center  of  mass  (i.e.,  the  lever  arm),  and  At  is  the  time  step. 
Since  the  location  of  th^  center  of  mass  is  not  known,  rj  is  not  known. 
The  displacement  of  jet  #j  from  the  center  of  mass  can  be  expressed 
as 


(5.3) 


where  rj'  is  the  known  location  on  the  vehicle  of  jet  #j,  and  Tcm  is  the 
location  of  the  center  of  mass  which  is  to  be  estimated. 

The  spacecraft  used  for  this  example  aerocapture  mission  is 
assumed  to  have  a  maximum  allowable  roll  rate  of  20  deg/sec.  With 
a  roll  rate  near  this  limit,  the  rate  change  due  to  the  torque  from  the 
dynamic  coupling  (i.e.,  cox  1(0)  is  significant  and  must  be  considered. 
Ignoring  this  term  for  now  in  order  to  put  Equation  5.2  in  a  form 
more  suitable  for  filter  design,  an  estimate  of  this  dynamic  coupling 
term  based  on  rate  measurements  and  the  current  inertia  matrix 
estimate  will  later  be  subtracted  from  the  difference  (i.e.,  the 
residual)  between  the  measured  and  predicted  rate  change. 

Without  dynamic  coupling.  Equation  5.2  can  be  rewritten  to 
give  the  rate  changes  experienced  by  a  spacecraft  due  to  jet  firings 
as 


96 


A©  =  r^a  -  r^(rcm  X  x)  +  dAt 


(5.4) 


where 


^  =  |XTjjAt  (5  5) 

is  defined  as  the  linear  impulse  due  to  all  jet  forces  and 

a  =  *'j  ^  (5.6) 

is  the  angular  impulse  of  the  vehicle  about  the  defined  vehicle  origin 
due  to  jet  firings.  The  disturbance  term  in  Equation  5.4  consists  of 
the  rate  changes  due  to  the  aerodynamic  forces,  the  unmodelled 
dynamics,  and  the  difference  between  the  actual  dynamic  coupling 
and  the  estimated  dynamic  coupling.  The  spacecraft  is  assumed  to 
operate  at  altitudes  high  enough  that  the  magnitudes  of  the 
aerodynamic  torques  are  small  compared  to  the  torques  applied  by 
the  RCS  jets.  Since  the  unmodelled  dynamics  and  the  difference 
between  the  actual  and  estimated  dynamic  coupling  terms  are  also 
considered  small,  the  last  term  of  Equation  5.4  can  be  considered  a 
noise  signal 


V  =  dAt  (5.7) 

With  known  jet  positions  and  thrust  levels,  the  only  unknowns 
in  Equation  5.4  other  than  the  noise  are  the  vehicle's  location  of 
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center  of  mass,  rcm»  and  the  inverse  inertia  matrix.  These  are  the 
mass  properties  required  by  the  actuator  selection  algorithm  to 
calculate  the  activity  vectors  for  each  jet.  With  a  noisy  system  like 
that  expressed  by  Equation  5.4,  an  estimator  such  as  a  Kalman  filter 
is  required  to  generate  accurate  parameter  estimates  in  the  presence 
of  noise. 

The  second-order  filter  designed  in  references  [23,24]  was 
derived  in  order  to  estimate  these  mass  properties  in  the  presence  of 
noise.  The  filter  contains  the  system  model  represented  by  Equation 
5.2  which  it  uses  to  predict  the  output  of  rate  gyros,  A^.  The 
predicted  output  is  then  used  in  comparison  with  the  actual 
measured  output,  and  revisions  to  the  mass  property  estimates  are 
made  accordingly. 

Since  the  filter  compares  measured  angular  rate  changes  to 
predicted  angular  rate  changes  at  discrete  times,  the  measured 
angular  rate  changes  must  first  be  expressed  as  a  discretized 
function  of  the  mass  properties  to  be  estimated.  A  state  vector  y 
whose  elements  are  the  mass  properties  to  be  estimated  is  defined  as 


y  =  In*  l22y  I33.  I23»  Ii3»  Ii2.  Temp  f, 


IT 


cm2»  TemaJ 


(5.8) 


where  is  the  element  in  the  i^h  row  and  column  of  the  inverse 
inertia  matrix,  and  remi  is  the  i*h  element  of  Fcm-  Rearranging 
Equation  5.4,  a  change  in  the  vehicle’s  angular  rates  can  be  expressed 
in  terms  of  y  as 


Ao)  =  h(y )  +  V 


(5.9) 
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where 


h(y)  =  Ay  +  y 

2  i=  1 

In  Equation  5.10,  A  is  a  3  x  9  matrix  satisfying 


Ay  =  r^ct  (5.11) 

while  each  Ai  is  a  symmetric  9x9  matrix  satisfying 

^y%y  =  [l‘MrcmXX)]i  (5.12) 

for  each  vector  component  i  =  1,  2,  3,  and  the  are  the  natural  basis 
vectors 


■  1  ■ 

O' 

o' 

5i  = 

0 

,  ^2  = 

1 

,  ^3  = 

0 

.0. 

.0. 

.  1 . 

(5.13) 


used  to  create  a  vector  out  of  the  scalar  terms,  y^Ajy.  The  matrices 
A  and  Ai  are  found  to  be: 


ai  0  0 

0  03  02 

0  0  0  ‘ 

A  = 

0  02  0 

03  0  oi 

0  0  0 

1 - 

o 

o 

p 

02  01  0 

0  0  0. 

(5.14) 
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(5.15a) 


Ai 


0  0  0 

0  0  0 

0  0  0 

0  0  0 

0.0  0 
0  0  0 

0  0  0 

-X3  0  0 

X2  0  0 


0  0  0  0  -X3  X2 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  ->.2  >.1  0 

0  0  0  A.3  0  -Xi 

0  -X2  X3  0  0  0 

0  Xi  0  0  0  0 

0  0  -Xi  0  0  0 


A2  = 


0 

0 

0 

0 

0 

0 

0 

0 

0 


0 

0 

0 

0 

0 

0 

^3 

0 

-^1 


0  0  0  0  0  0  0 

0  0  0  0  X3  0  -A.1 

0  0  0  0  0  0  0 

0  0  0  0  -X2  0 

0  0  0  0  0  0  0 

0  0  0  0  0  -A,3  X2 

0  -X2  0  0  0  0  0 

0  0  -A,3  0  0  0 

0  0  0  X2  0  0  0 


(5.15b) 
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(5.15c) 
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0 
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0 
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0 
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0 

0 

0 

0 

0 

0 

0 

-X.1 

-A,2 

0 

0 

0 

0 

If  only  the  elements  of  the  inverse  inertia  matrix  were 
required,  the  linear  portion  of  Equation  5.10  would  suffice;  however, 
the  location  of  the  center  of  mass  is  also  needed.  Since  the  third  3  x 
3  block  partition  of  A  contains  only  zeros,  the  elements  of  Fcm  do  not 
contribute  to  the  linear  portion  of  Equation  5.10.  The  complete 
nonlinear  model  represented  by  Equation  5.10  must  therefore  be 
used  in  order  to  estimate  Fcm. 

Neglecting  the  effects  of  dynamic  coupling,  the  difference  in  the 
discrete  measurements  of  the  angular  rates  from  time  k-1  to  time  k 
is  assumed  to  satisfy  the  stochastic  equation 

A(0n,  =  h{yk)  +  Vk  (5.16) 

where  Acum  is  the  measured  output  of  the  gyros,  h{yk)  is  from 
Equation  5.10  at  time  k,  and  Vk,  consisting  of  both  the  disturbance 
noise  V  and  the  sensor  noise,  is  assumed  to  be  a  3-dimensional  white 
gaussian  noise  signal  with  zero  mean  and  independent  of  Vk  and  all 
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forces  and  torques.  This  assumed  gaussian  distribution  for  the 
combined  noise  signal  is  considered  reasonable  since  the  sensor  noise 
is  assumed  to  dominate  and  common  gyros,  such  as  the  DRIRU  II 
[24],  exhibit  gaussian  noise  characteristics. 

Successive  updates  to  the  mass  property  estimates  are  then 
generated  by  the  estimator  from 


yk  =  yk-i  + 


(5.17) 


where  Kk  is  a  gain  matrix  calculated  via  equations  2-119  through  2- 
121  of  reference  [24]  at  each  step  and  is  a  function  of  the  error 
covariance,  the  jets  being  fired,  and  the  noise  characteristics  of  the 
sensors;  AtOm  is  the  measured  difference  in  angular  rates  from  time 
k-1  to  time  k;  Ao^  is  the  predicted  angular  rate  change  determined 
from  the  rigid  spacecraft  model  of  Equation  5.2;  and  I  (wm  X  iwm)  is 
the  estimated  rate  change  due  to  the  dynamic  coupling  effects. 

5.3  Input  Selection  During  Dedicated  Mass  Property 
Identification 

The  filter  design  described  above  generates  mass  property 
estimates  based  on  the  comparison  of  predicted  and  measured 
angular  rate  changes  due  to  RCS  jet  firings.  In  order  to  rapidly 
converge  to  accurate  estimates,  however,  these  jet  firings  cannot  be 
random.  For  example,  little  can  be  learned  about  a  vehicle's  roll 
inertia  by  firing  only  pitch  or  yaw  jets.  A  jet  selection  strategy 
which  chooses  the  jet  firings  that  are  most  likely  to  maximize  the 
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increase  in  knowledge  of  the  estimated  parameters  was  originally 
developed  in  reference  [35].  This  strategy  is  employed  in  the 
candidate  design  during  the  periods  dedicated  to  mass  property 
estimation. 

The  goal  of  this  input  selection  strategy  is  to  develop  an 
analytical  method  to  determine  the  maneuver  (i.e.,  jet  firings)  that 
will  maximize  the  increase  in  knowledge  of  the  desired  mass 
properties.  To  do  so,  the  error  covariance  matrix  propagation  across 
a  measurement  is  examined.  The  change  in  the  covariance  matrix 
due  to  a  measurement  at  time  k  is  given  by 

APk  =  Kk  Hk  Pk(-)  (5.18) 

where  Kk  is  the  gain  matrix  at  time  k,  Hk  is  the  output  matrix  at  time 
k,  and  Pk(-)  is  the  error  covariance  matrix  prior  to  the  measurement. 
A  detailed  explanation  of  the  terms  in  Equation  5.8  can  be  found  in 
reference  [24]. 

Equation  5.18  depends  only  on  the  old  state  estimate,  the  old 
error  covariance  matrix,  the  noise  covariance  matrix,  and  the  current 
jet  thrusts  and  positions  (i.e.,  the  applied  torques).  Since  all  of  these 
are  known  or  have  estimates  (whether  accurate  or  not).  Equation 
5.18  can  be  evaluated  prior  to  each  step  to  determine  the  increase  in 
accuracy  of  the  mass  property  estimates  due  to  the  firing  of  a 
particular  jet.  A  decrease  in  the  trace  of  P  represents  an  increase  in 
accuracy  of  the  mass  property  estimates.  At  each  step,  the  jet  which 
maximizes  the  increase  in  accuracy  (i.e.,  maximizes  the  decrease  in 
the  trace  of  P)  is  then  chosen  to  be  fired.  The  chosen  jet  firings 
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result  in  a  sequence  which  maximizes  the  likelihood  of  rapid 
convergence  of  the  mass  property  estimates,  thus  minimizing  fuel 
consumption  during  the  estimation  process. 
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Chapter  6 


System  Integration 

6.1  Introduction 

Chapters  3,  4,  and  5  presented  the  development  of  the  three 
main  elements  of  the  candidate  aerocapture  attitude  control  system: 
the  robust  control  law,  the  jet  selection  algorithm,  and  the  mass 
property  identification  algorithm.  When  combining  these  elements 
into  the  final  attitude  control  system  structure  first  illustrated  in 
Figure  2.1,  several  integration  issues  must  be  addressed.  Such  issues 
include  the  system  cycling  rates,  the  interaction  of  the  main 
elements,  and  the  timing  of  the  mass  property  updates.  The  overall 
system  cycling  rate,  and  the  cycling  rates  of  each  of  the  main 
elements,  are  discussed  in  Section  6.2.  Section  6.3  describes  the 
interaction  of  the  control  law  and  the  jet  selection  algorithm  along 
with  the  control  law  "tuning”  required  to  make  these  elements  work 
together.  Section  6.4  describes  how  the  mass  property  identification 
algorithm  is  employed  to  detect  mass  property  changes  and  provide 
timely  updates  to  the  mass  property  estimates.  A  brief  summary  of 
the  overall  operation  of  the  candidate  aerocapture  attitude  control 
system  is  provided  in  Section  6.5.  Section  6.6  outlines  how  the 
candidate  attitude  control  system  can  be  applied  to  attitude  control 
for  other  aerocapture  vehicles 
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6.2  System  Cycling  Rates 

Since  the  candidate  control  system  was  designed  for  a  mission 
and  vehicle  similar  to  the  AFE,  the  decision  as  to  what  rates  to  cycle 
the  overall  attitude  control  system  and  its  elements  was  based  on 

cycling  rates  of  the  existing  AFE  DAP.  Cycling  rates  of  the  control 
system  for  use  in  other  aerocapture  vehicles  may  vary  and  will 

usually  be  dictated  by  system  specifications. 

The  AFE  guidance  law  generates  bank  angle  commands  once 
per  second  (i.e.,  1  hz).  The  existing  AFE  DAP  operates  at  a  25  hz  rate 

corresponding  to  the  minimum  RCS  jet  firing  time  of  40  msec.  The 

AFE  DAP’s  proportional  control  law  and  its  table  look-up  jet  selection 
algorithm  are  both  cycled  at  this  25  hz  rate.  Each  time  the  control 
law  is  cycled  it  generates  angular  acceleration  commands  which  are 
used  in  the  table  look-up  jet  selection  algorithm  to  determine  which 
jets  to  fire. 

For  the  candidate  control  system,  the  guidance  is  also  assumed 
to  provide  bank  angle  commands  at  a  1  hz  rate.  Like  the  AFE  DAP, 
the  sliding  mode  control  law  is  cycled  at  a  25  hz  rate.  The  linear 
programming  jet  selection  algorithm,  however,  is  cycled  at  2.5  hz. 
The  reason  for  the  lower  rate  is  due  to  the  fact  that  the  solution  of 
each  linear  program  consists  of  the  jet  duty  cycles  required  to 
provide  the  angular  accelerations  commanded  by  the  control  law.  In 
implementation,  the  400  msec  between  each  solution  to  the  linear 
program  is  divided  into  ten  periods  of  40  msec  in  length.  Jets  are 
fired  for  the  fraction  of  these  ten  periods  corresponding  to  their 
associated  decision  variables  from  the  last  solution  of  the  linear 
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program.  For  example,  if  the  last  solution  of  the  linear  programming 
problem  results  in  a  decision  variable  (i.e.,  commanded  duty  cycle)  of 
xj=0.41  for  jet  #j,  jet  #j  will  be  fired  during  four  of  the  ten  40  msec 
time  periods  between  solutions  to  the  linear  program. 

To  more  evenly  implement  the  commanded  duty  cycles  over 
the  ten  periods  between  solutions  to  the  linear  program,  a  running 
ratio  is  calculated  for  each  selected  jet  to  determine  if  the  jet  should 
be  fired  during  the  present  period  [21].  The  ratio  is 


Ni 


0.5+  X  Jj(k) 

k=  1 
Nt 


(6.1) 


where  Nt  is  the  number  of  periods  that  have  elapsed  since  the  last  jet 
selection,  and  Jj  is  defined  by 


1  if  jet  #j  was  fired  during  time  cycle  #k 
0  otherwise 


(6.2) 


A  selected  jet  is  fired  for  the  current  40  msec  time  period  if  the 
following  condition  is  true: 


^jk<Xj  (6.3) 

Table  6.1  shows  an  example  of  how  the  condition  6.3  is  used  to 
evenly  spread  jet  firings  throughout  the  ten  time  periods  between 
solutions  of  the  jet  selection  linear  program.  The  entries  in  the 
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second  column  of  Table  6.1  represent  the  calculation  of  the  ratio 
(Equation  6.1)  for  a  jet  with  a  selected  decision  variable  of  xj  =  0.41. 
The  third  column  contains  the  decisions  of  whether  or  not  to  fire  the 
jet  based  on  Equation  6.3. 


time  step 

fire  jet 
(Ri  <  xi)? 

1 

0.50 

no 

2 

0.25 

yes 

3 

0.50 

no 

4 

0.38 

yes 

5 

0.50 

no 

6 

0.42 

no 

7 

0.36 

yes 

8 

0.44 

no 

9 

0.39 

yes 

10 

0.45 

no 

Table  6.1.  Example  of  Firing  of  Jets  Between  Solutions  of 
Jet  Selection  Linear  Program  (xj  =  0.41) 

Jet  #j  is  fired  for  four  of  the  ten  40  msec  periods  between  solutions 
of  the  linear  programming  problem. 

Like  the  control  law,  the  mass  property  identification  algorithm 
is  also  cycled  at  25  hz.  Identification  is  performed  on  every  jet 
firing,  but  updated  mass  property  estimates  are  provided  to  the 
control  law  and  jet  selection  algorithm  only  when  needed.  Section 
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6.4  outlines  the  scheme  for  determining  when  to  provide  these 
updates. 

6.3  Control  Law/Jet  Selection  Interaction 

The  ability  of  the  jet  selection  algorithm  to  provide  the  angular 
accelerations  commanded  by  the  control  law  is  limited  by  the  control 
authority  of  the  jets.  Since  the  control  law  generates  the  inputs  to 
the  jet  selection  algorithm  without  knowledge  of  the  limitations  of 
the  actuators,  it  is  possible  for  the  control  law  to  command  angular 
accelerations  that  the  jet  selection  algorithm  cannot  provide.  From 
the  equations  for  the  sliding  surface  (Equation  3.32)  and  the  sliding 
condition  (Equation  3.31),  it  is  obvious  the  control  law  parameters 
and  r\i  determine  how  quickly  the  control  law  attempts  to  drive  the 
vehicle  to  its  desired  state.  If  these  parameters  are  chosen  so  that 
they  dictate  a  response  quicker  than  the  vehicle  is  capable  of 
providing,  the  control  law  may  consistently  command  accelerations 
greater  than  the  jet  selection  algorithm  can  provide.  When  this 
happens,  the  actual  accelerations  of  the  vehicle  will  lag  the 
commanded  accelerations  and  the  vehicle  will  have  difficulty  in 
tracking  the  desired  trajectory. 

This  problem  is  demonstrated  in  a  simulation  of  the  combined 
control  law/jet  selection  algorithm.  In  Chapter  3,  each  was  chosen 
from  the  rule  of  thumb 
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to  be  5.  Since  the  only  requirement  on  each  qi  is  that  they  be 
positive  constants,  let  each  q,  be  chosen  as  1.  Figure  6.1  shows  the 
bank  angle  tracking  response  of  the  nominal  AFE  vehicle  using  the 
new  control  law  and  jet  selection  algorithm  during  a  simple  120 
second  simulation.  The  simulation  is  started  with  the  vehicle  at  0 
deg  bank  and  trim  a  and  p.  The  commanded  bank  angle  is  initially 
45  deg  and  varies  sinusoidally  during  the  simulation.  A  commanded 
bank  reversal  occurs  60  sec  into  the  simulation. 
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Figure  6.1.  Bank  Angle  Tracking 
(^i=  5.0,  q,'  =  1.0) 


With  Xi  =  5,  T|i  =  1  and  an  initial  45  deg  offset.  Equation  3.12 
predicts  that  the  sliding  surface  should  be  reached  (resulting  in 
tracking  to  within  2  deg)  in  a  time  less  than 

I  Si(x(0),  O) I  _  I  Xj [({>(0)  -  (})d(0)]  I  _  ^  -  3  9  ^ec 

1  ■ 

Figure  6.1  shows  that  the  vehicle  attempts  to  reach  the  desired  bank 
angle  trajectory  rapidly,  but  is  not  able  to  do  so  as  quickly  as 
dictated  by  Equation  6.5  due  to  the  limited  authority  of  the  RCS  jets. 
Also,  because  of  this  limited  control  authority,  the  vehicle  is  not  able 
to  reverse  its  velocity  as  quickly  as  the  control  law  dictates,  and  the 
vehicle  overshoots  the  desired  trajectory.  When  correcting  the 
tracking  error  due  to  the  overshoot,  the  vehicle  again  overshoots  the 
desired  trajectory  in  the  other  direction.  Though  the  state  initially 
reaches  the  sliding  surface  rapidly,  it  does  not  settle  to  the  desired 
trajectory  for  another  30  seconds. 

Since  the  jets  cannot  provide  the  needed  control  authority  for 
the  vehicle  to  react  as  fast  as  required  by  Equation  6.5,  the  control 
law  must  be  "tuned"  by  altering  the  Xj's  and  "ni’s  so  that  the  control 
law  does  not  consistently  command  angular  accelerations  the  jet 
selection  process  cannot  provide.  Using  the  AFE  DAP  in  the  same 
simulation  as  above,  the  AFE  vehicle  is  able  to  control  to  the  desired 
trajectory  with  an  initial  45  deg  offset  in  approximately  10  sec. 
Using  10  sec  as  a  goal,  the  simulation  above  was  performed 
numerous  times  with  various  combinations  of  X[  and  r\\.  Using  this 


trial -and-error  approach,  a  combination  of  Xi  =  1.0  and  r\i  =0.1  is 
found  to  give  the  desired  response  with  a  minimum  of  jet  firings. 
With  this  combination,  the  sliding  surface  should  be  reached 
(resulting  in  tracking  to  within  2  deg)  in  a  time  less  than 

Tli  0.1 

Figure  6.2  shows  that,  with  the  new  parameters,  the  actual  bank 
angle  trajectory  reaches  to  within  2  deg  of  the  desired  trajectory 
within  8.0  sec  and  that  overshoot  is  eliminated.  With  A.i  =  1.0  and  ti, 
=  0.1,  the  control  law  is  properly  tuned. 

This  control  law  tuning  may  also  be  required  when  using  this 
candidate  attitude  control  system  design  for  other  aerocapture 
vehicles.  The  actuator  control  authority  of  the  particular  vehicle  will 
dictate  the  required  tuning. 
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Figure  6.2.  Bank  Angle  Tracking 
(X|=  1.0,  m  =  0.1) 


6.4  Mass  Property  Updates 

From  Chapter  4,  the  linear  programming  jet  selection  algorithm 
characterizes  each  jet  available  for  selection  by  an  activity  vector 
(i.e..  Equation  4,3)  which  represents  the  effect  of  firing  the  jet. 
Accuracy  of  the  activity  vectors  is  limited  by  errors  in  the  estimates 
of  the  vehicle  inertia  matrix  and  center  of  mass  location.  The  mass 
property  identification  algorithm  must  provide  accurate  initial 
estimates  of  these  mass  properties  and  must  update  the  estimates  if 
the  mass  properties  change  significantly  during  the  mission.  Failure 


to  provide  accurate,  timely  estimates  may  result  in  the  aerocapture 
attitude  control  system  attempting  to  control  a  vehicle  with  mass 
properties  outside  the  range  for  which  it  is  robust.  When  this 
happens,  the  vehicle  may  become  unstable.  This  section  describes 
how  the  mass  property  identification  algorithm  is  employed  to 
provide  both  the  initial  mass  property  estimates,  and  any  required 
updates. 

6.4.1  Initial  Estimation 

The  method  of  employing  the  mass  property  identification 
algorithm  for  the  initial  estimates  is  straight  forward.  The  estimates 
of  the  mass  properties  to  be  identified  are  initially  set  to  zero  in  the 
extended  Kalman  filter.  For  a  period  of  0.5  sec  beginning  at  the  point 
of  atmosphere  entry  at  an  altitude  of  400,000  feet,  the  control  law 
and  the  jet  selection  algorithm  are  bypassed  and  vehicle  motions  are 
dedicated  strictly  to  mass  property  identification.  Operating  at  25  hz, 
a  sequence  of  12  jet  firings  is  generated  by  the  selection  strategy 
outlined  in  Section  5.3.  Comparing  the  predicted  angular  rate 
changes  to  the  measured  rate  changes  due  to  these  12  firings,  the 

estimator  is  able  to  converge  to  accurate  mass  property  estimates  by 
the  end  of  the  dedicated  identification  period.  At  the  end  of  the 

initial  identification  period,  the  mass  property  estimates  are  supplied 
to  both  the  control  law  and  the  jet  selection  algorithm.  Activity 
vectors  for  each  jet  are  calculated  in  the  jet  selection  algorithm,  and 

the  control  law  and  jet  selection  algorithm  begin  to  implement  the 

bank  angle  commands  from  the  guidance  law. 
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6.4.2  Detecting  Mass  Property  Changes 

The  control  law  and  jet  selection  algorithm  alone  represent  a 
robust  control  system.  Small  mass  property  changes  can  be  tolerated 
without  adversely  affecting  the  control  system  performance.  Since 
the  system  may  become  unstable  with  large  mass  property  changes, 
the  control  system  must  somehow  detect  significant  changes  in  the 
mass  properties.  After  detection,  mass  property  identification  can  be 
performed  to  maintain  stability. 

What  determines  a  "significant”  mass  property  change?  Any 
mass  property  change  that  puts  the  mass  properties  outside  the 
range  for  which  the  control  system  is  robust  causing  the  vehicle  to 
become  unstable  is  certainly  significant  and  must  be  detected. 
Additionally,  for  mass  properties  that  change  slowly,  it  is  necessary 
to  detect  these  changes  before  the  vehicle  becomes  unstable. 

During  the  normal  operation  of  the  attitude  control  system,  the 
mass  property  identification  algorithm  continually  monitors  the 
effects  of  jet  firings  on  the  vehicle  dynamics  in  order  to  detect 
significant  mass  property  changes.  For  each  jet  firing,  the  estimator 
predicts  the  resulting  angular  rate  changes  using  the  last  mass 
property  estimates  supplied  to  the  control  law  and  jet  selection 
algorithm.  When  the  magnitude  of  the  residual  (i.e.,  difference) 
between  the  estimated  rate  changes  and  the  rate  changes  measured 
by  the  spacecraft  gyros  becomes  large  enough  to  indicate  a 
significant  mass  property  change,  a  flag  is  set  signaling  that  updated 
mass  properties  estimates  are  required. 


The  residual  value  used  to  trigger  the  flag  is  found  by 
assuming  the  control  law/jet  selection  algorithm  combination  can 
tolerate  up  to  a  40  inch  displacement  in  the  center  of  mass  location 
with  the  associated  change  in  the  inertia  matrix  calculated  as 
described  in  Appendix  C  (this  assumption  will  be  shown  to  be 
reasonable  from  the  robustness  test  results  presented  in  Chapter  7). 
With  a  40  inch  displacement,  the  center  of  mass  is  located  at 

fcm  — 

and  the  inertia  matrix  is  calculated  to  be 

■  4069.79  -  915.79  -  974.39 

I  =  -  915.79  3280.79  -  889.79 

.  -  974.39  -  889.79  3000.19 


slugs  fr 


(6.8) 


106.69 
23.06 
L  22.48  J 


inches 


(6.7) 


The  residual  between  the  measured  and  predicted  angular  rate 
change  due  to  the  firing  of  RCS  jets  is  defined  as 

R  =  A(0ni-A©p  (6.9) 

where 

R  =  the  residual  vector 

AOm  =  measured  angular  rate  changes 

AcOp  =  predicted  angular  rate  changes 


Using  Equation  5.2  for  the  expression  of  the  rate  changes,  the 
residual  due  to  the  firing  of  a  single  jet  #j  can  be  written  as 


R  =  r^ 


W  X  1(0  +  (Fj  -  Fcm)  X  Tj  +  d 
W  X  1(0  +  (fj  -  fcm)  X  Tj  +  d 


At 

At 


(6.10) 


The  nominal  AFE  vehicle  mass  properties  [36]  are  listed  in 
Table  6.2  with  the  nominal  inertia  matrix  calculated  as  in  Appendix  B 
relative  to  the  body  axes  with  origin  at  the  vehicle  center  of  mass. 
The  position  and  thrust  vector  for  each  of  the  AFE's  16  RCS  jets  [31] 
are  listed  in  Table  6.3  with  the  coordinate  system  used  to  specify  the 
location  of  the  center  of  mass  illustrated  in  Figure  6.3. 


inertia 

r  2248.2  -5.0  -63.6  1 

1  =  -5.0  1459.2  21.0 

L  -63.6  21.0  1178.6. 

slugs  •  ft^ 

center  of  mass 

r  83.61 

fcm  =  -0.03 

L-O.61J 

inches 

mass 

mass  =  123.0 

slugs 

Table  6.2,  Nominal  Vehicle  Mass  Properties 


Jet 

# 

Location  (inches) 

Thrust  Vector 

(Ibf)  1 

X 

y 

z 

X 

y 

z 

1 

60.69 

-6.5 

-45.0 

0 

+30 

0 

2 

60.69 

+6.5 

+45.0 

0 

-30 

0 

3 

60.69 

-45.0 

-6.5 

0 

0 

+30 

■1 

60.69 

+45.0 

+6.5 

0 

0 

-30 

5 

60.69 

+3.5 

+45.0 

+30 

0 

0 

6 

60.69 

-3.5 

-45.0 

+30 

0 

0 

■a 

60.69 

-45.0 

-3.5 

+30 

0 

0 

8 

60.69 

+45.0 

+3.5 

+30 

0 

0 

9 

60.69 

0 

+41.5 

+  125 

0 

0 

10 

60.69 

0 

+48.5 

+  125 

0 

0 

1  1 

60.69 

0 

-41.5 

+  125 

0 

0 

12 

60.69 

0 

-48.5 

+  125 

0 

0 

13 

60.69 

-41.5 

0 

+  125 

0 

0 

14 

60.69 

-48.5 

0 

+  125 

0 

0 

1  5 

60.69 

+41.5 

0 

+  125 

0 

0 

16 

60.69 

+48.5 

0 

+  125 

0 

0 

Table  6.3.  Vehicle  Jet  Properties 
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Figure  6.3.  Body  Referenced  Coordinate  System 


The  typical  magnitude  of  the  residual  vector  associated  with 
the  firing  of  a  single  jet  is  used  as  the  trigger  value.  Since  the 
majority  of  the  jet  firings  commanded  by  the  jet  selection  algorithm 
are  roll  jets,  jet  #4  is  chosen  to  determine  the  typical  residual  vector 
magnitude.  Assuming  there  are  no  disturbances  or  angular  rates 
(i.e.,  d  =  (0  =  0),  the  estimated  mass  properties  are  the  nominal  values 
listed  in  Table  6.2,  the  actual  mass  properties  are  those  for  a  40  inch 
displacement  given  by  Equations  6.7  and  6.8,  and  the  thrust  and 
position  of  jet  #4  are  as  listed  in  Table  6.3,  the  residual  vector  due  to 
a  single  firing  of  jet  #4  is  calculated  from  Equation  6.10  to  be 


(6.11) 


R  = 


0.7768 

-0.4409 

.-0.9179. 


X  10'^  tads, 
sec 


The  magnitude  of  this  residual  vector  is  1.3  x  10" 3  rads/scc.  To 
ensure  a  mass  property  change  is  detected  before  this  large  a  change 
occurs,  the  detection  trigger  value  is  chosen  to  be  1.0  x  10*3  rads/sec. 
Using  this  approach,  a  trigger  value  based  on  the  jet  and  mass 
properties  of  other  vehicles  can  be  calculated  for  other  aerocapture 
missions. 

As  jets  are  fired,  the  residual  vectors  are  calculated  and  the 
magnitudes  are  stored.  To  prevent  triggering  on  a  single  residual 
measurement  which  may  jump  due  to  noise,  the  average  of  the  last 
three  residual  magnitudes  is  calculated.  The  flag  indicating  that  a 
significant  mass  property  change  has  occurred  is  then  set  when  this 
average  becomes  larger  than  the  1.0  x  10*3  rads/sec  trigger  value. 


6.4.3  Updating  Mass  Property  Estimates  After  a  Change  is 
Detected 

When  the  detection  flag  is  set,  the  extended  Kalman  filter  is 
reinitialized,  and  the  identification  of  the  new  mass  properties  is 
begun.  New  mass  property  estimates  are  not  supplied  to  the  control 
law  and  jet  selection  algorithm  until  the  end  of  the  new  identification 
period.  Unlike  the  initial  dedicated  identification  period,  the  control 
law  and  jet  selection  algorithm  are  not  bypassed  after  a  change  is 
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detected.  Instead,  in  an  effort  to  reduce  the  impact  of  the 
identification  process  on  the  attitude  control  task,  the  identification 
algorithm  attempts  to  estimate  the  new  mass  properties  based  on  the 
jet  firings  commanded  by  the  jet  selection  algorithm. 

Jet  firings  in  addition  to  those  commanded  by  the  jet  selection 
algorithm  are  required.  The  overall  attitude  control  task  is  to  track 
the  bank  angle  trajectory  commanded  by  the  guidance  law.  Bank 
maneuvers  are  performed  primarily  by  a  combination  of  roll  and 
yaw  maneuvers;  pitch  jets  are  rarely  commanded  and  little 
information  about  the  inertia  about  the  pitch  axis  can  be  gained.  To 
alleviate  this  problem,  the  selection  scheme  used  during  the  initial 
identification  period  is  invoked  during  one  out  of  every  five  cycles 
(i.e.,  5  hz  instead  of  25  hz).  The  one  jet  chosen  by  this  selection 
scheme  is  then  fired  along  with  any  jets  commanded  by  the  jet 
selection  algorithm. 

The  new  identification  period  continues  until  the  following  two 
conditions  are  met.  The  first  condition  is  that  the  filter  estimates 
must  improve  the  average  residual  magnitude  by  at  least  an  order  of 
magnitude  (i.e.,  becomes  less  than  1.0  x  lO*'^  rads/sec).  The  second 
condition  is  that  all  the  diagonal  elements  of  the  error  covariance 
matrix  P  reduce  to  within  an  order  of  magnitude  of  what  they  were 
before  the  detected  change.  The  first  condition  ensures  a  minimum 
acceptable  accuracy  of  the  estimates,  while  the  second  condition 
ensures  the  identification  period  is  not  ended  while  the  filter  is  still 
converging  rapidly  to  more  accurate  estimates. 
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Like  any  Kalman  filter  design,  after  several  updates  to  the 
estimated  state  vector  (i.e..  Equation  5.17)  the  filter  becomes 
optimistic  and  the  information  from  later  jet  firings  is  weighed  less 
than  earlier  firings.  To  prevent  the  filter  from  becoming  optimistic, 
the  elements  of  the  error  covariance  matrix  are  set  artificially  high  at 
the  beginning  of  a  new  estimation  period.  As  the  information  from 
more  and  more  jet  firings  is  extracted  and  the  estimator  converges  to 
accurate  estimates  of  the  new  mass  properties,  the  dements  of  the 
error  covariance  matrix  become  smaller.  By  the  time  the  elements  of 
P  become  small  enough  to  meet  the  second  condition,  little 

information  can  be  extracted  from  additional  jet  firings  and  the  filter 

should  have  converged  to  accurate  estimates  of  the  new  mass 
properties. 

Once  the  identification  period  ends,  the  new  mass  property 
estimates  are  supplied  to  the  control  law  and  the  jet  selection 
algorithm.  New  activity  vectors  for  each  jet  are  calculated  in  the  jet 
selection  algorithm,  and  the  mass  property  identification  algorithm 
again  begins  monitoring  the  magnitude  of  the  residual  vector  to 
detect  any  further  mass  property  changes. 

6.5  Summary  of  Operation 

Figures  6.4,  6.5,  and  6.6  depict  the  operation  of  the  candidate 
control  system  during  the  three  main  modes  of  employment.  As 
Figure  6.4  illustrates,  the  control  law  and  jet  selection  algorithm  are 
bypassed  during  the  initial  dedicated  mass  property  identification 
algorithm.  During  this  first  0.5  sec  of  operation,  which  occurs  at  the 


point  of  atmosphere  entry,  the  motion  of  the  vehicle  is  dedicated  to 
mass  property  identification  with  the  input  selection  scheme  in  the 
mass  property  identification  algorithm  commanding  the  jet  firings 
needed  to  rapidly  converge  to  accurate  mass  property  estimates. 


Figure  6.4.  Dedicated  Mass  Property  Identification  Period 

(First  0.5  sec) 
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Figure  6.5.  Normal  Control  System  Operation 


After  the  initial  identification  period,  the  control  system  begins 
its  normal  operation  with  the  control  law  and  the  jet  selection 
algorithm  being  used  to  implement  the  bank  angle  trajectory 
commanded  by  the  guidance  law  (see  Figure  6.5).  The  control  law 


operates  at  25  hz,  the  linear  programming  jet  selection  problem  is 
solved  at  2.5  hz,  and  jets  are  commanded  on  or  off  at  25  hz  to 

implement  the  solutions  to  the  linear  program.  During  this  normal 

operation,  the  mass  property  identification  algorithm  monitors  the 
effects  of  jet  firings  on  the  rate  changes  experienced  by  the  vehicle 
and  compares  the  measured  rate  changes  to  its  predicted  (i.e., 
modelled)  rate  changes.  When  the  magnitude  of  the  residual  vector 

becomes  larger  than  the  trigger  value  of  1.0  x  10'3  rads/sec,  a  flag  is 

set  indicating  a  change  in  mass  properties  has  been  detected.  When 
this  flag  is  set,  normal  control  system  operation  ends  and  a  new  mass 
property  identification  period  begins. 

During  the  new  identification  period,  the  control  law  and  the  jet 
selection  algorithm  continue  their  task  of  implementing  the 
commanded  bank  angle  trajectory.  As  Figure  6.6  illustrates,  mass 
property  identification  is  attempted  based  on  the  jet  firings 
commanded  by  the  jet  selection  algorithm^  In  addition  to  these  jet 
firings,  the  input  selection  scheme  described  in  Section  5.3  is  invoked 
at  a  5  hz  rate  and  is  used  to  choose  a  jet  firing  to  help  speed  the 
convergence  of  the  identification  algorithm  to  accurate  mass 
property  estimates.  This  new  identification  period  continues  until 
the  two  conditions  described  in  Section  6.4.3  are  met. 

At  the  end  of  the  new  identification  period,  the  updated  inass 
property  estimates  are  supplied  to  the  control  law  and  the  jet 
selection  algorithm.  New  activity  vectors  based  on  the  new  mass 
property  estimates  are  calculated  in  the  jet  selection  algorithm,  and 
the  attitude  control  system  returns  to  its  normal  operation  with  the 
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mass  property  identification  algorithm  monitoring  the  magnitude  of 
'•  '  resi'^''al  vectors  for  detection  of  additional  mass  property 

Cl 


•onsad 

aeceteratlon 


MASS  PROP  ID 

(25  hz) 


estimated  velocity  angles 


Figure  6.6.  Operation  After  Detection  of  Mass  Property 

Change 
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6.6  Application  to  Other  Aerocapture  Missions 

The  candidate  aerocapture  control  system  structure  outlined  in 
Figures  6.4,  6.5,  and  6.6  was  designed  to  be  used  for  attitude  control 

of  any  aerocapture  vehicle  or  mission.  The  overall  structure  is  also 

transportable  to  other  missions  or  vehicles.  The  following  mission 
and  vehicle  dependent  issues,  however,  must  be  addressed: 

A.  Control  Law 

(1)  The  bounds  on  the  unknowns  in  the  system  model  (i.e.. 

Equation  3.43)  must  be  calculated  as  in  Appendix  B  based 
on  the  vehicle  mass  and  aerodynamic  properties,  along 

with  the  atmospheric  properties  expected  during  the 
mission. 

(2)  The  choices  of  Xi  and  qj  must  be  made  based  on  the 

system  cycling  rates  and  the  control  authority  of  the 

vehicle's  jets.  Control  law  tuning,  as  discussed  in  Section 
6.3,  may  be  required. 

B.  Linear  Programming  Jet  Selection 

(1)  Activity  vectors  for  each  jet  (i.e..  Equation  4.3)  depend  on 
the  vehicle  inertia  matrix  and  center  of  mass  location 
provided  by  the  mass  property  identification  algorithm. 
They  also  depend  on  the  location  and  thrust  vector  of 
each  jet.  These  jet  properties  must  be  provided  for  the 
vehicle  to  be  used. 

(2)  The  relative  cost  of  using  each  jet  must  be  specified 
based  on  their  fuel  consumption  rates. 

(3)  A  strategy,  such  as  described  in  Section  6.2,  for 
commanding  jets  on  and  off  to  implement  the 
commanded  duty  cycles  must  be  developed  based  on 
system  cycling  rates  and  minimum  jet  on  times. 
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C.  Mass  Property  Identification 

(1)  Measurement  noise  covariance  must  be  determined  as  in 
reference  [24]  based  on  the  noise  characteristics  of  the 
vehicle  sensors. 

(2)  The  magnitude  of  the  residual  vector  used  to  trigger  the 
flag  indicating  detection  of  a  mass  property  change  must 

■  be  calculated  based  on  the  vehicle  mass  and  jet 
properties  as  in  Section  6.4.2. 
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Chapter  7 


Simulation  Results 

7.1  Introduction 

The  candidate  aerocapture  attitude  control  system  has  been 
tested  in  a  computer  simulation  of  an  Earth  aerocapture.  Section  7.2 
briefly  describes  the  simulation  developed  for  this  effort.  Results  of 
testing  to  demonstrate  the  increased  robustness  of  the  new  design 
over  the  existing  AFE  DAP  design  are  included  in  Section  7.3.  Section 
7.4  presents  results  of  tests  demonstrating  the  ability  of  the 
candidate  system  to  identify  and  adapt  to  changing  mass  properties. 
The  tests  described  in  Section  7.4  are  performed  without  modelled 
aerodynamic  forces  and  torques.  Section  7.5  describes  the  effects 
aerodynamics  have  on  these  tests. 

7.2  Simulation  Description 

A  6-DOF  simulation  of  an  earth  aerocapture,  written  in 
FORTRAN,  was  developed  for  testing  purposes  and  is  hosted  on  a 
Macintosh  Ilex  personal  computer.  It  simulates  the  actual  AFE 
mission  for  a  spacecraft  similar  to  the  AFE  vehicle,  but  with  a  wider 
range  of  possible  mass  properties.  For  the  purposes  of  this  thesis, 
the  AFE  guidance  has  not  been  included  in  the  simulation.  Instead, 
simple  bank  angle  commands  mimicking  the  typical  characteristics  of 
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a  trajectory  commanded  by  the  AFE  guidance  are  generated  and 
supplied  to  the  attitude  control  system  during  simulation. 


7.2.1  Initial  Conditions 

For  each  test,  the  simulated  aerocapture  maneuver  begins  with 
the  spacecraft  entering  Earth's  atmosphere  with  the  following  entry 
conditions  matching  those  projected  for  the  AFE  mission  [lOJ: 


altitude  =  400,000  ft 
latitude  =  0  deg 
longitude  =  0  deg 

velocity  magnitude  =  33,820  ft/sec 
VI  (North)  =  0 
V2  (East)  =  33,716  ft/sec 
V3  (toward  center  of 

Earth)  =  2,653  ft/sec 
flight  path  angle  =  4,5  deg  below 

local  horizontal 


At  entry,  the  vehicle  is  in  the  following  attitude: 

bank  angle  i^)  =  0  deg  (7.2) 

angle-of-attack  (a)  =  17.0  deg 
sideslip  angle  (P)  =  0  deg 

with  <[)  =  0  deg  indicating  full  lift  up  (i.e.,  lift  in  the  vertical  plane  only 
and  directed  away  from  the  center  of  the  Earth). 
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7.2.2  Vehicle  Properties 

For  comparison  purposes,  simulations  using  the  nominal  mass 
properties  listed  in  Table  6.2  are  run  with  both  the  AFE  DAP  and  the 
candidate  attitude  control  system.  Additional  simulations  are  run 
with  the  vehicle  mass  properties  varying  from  these  nominal  values. 
The  methods  used  to  vary  the  mass  properties  are  described  in 
Sections  7.3  and  7.4  below  for  each  test. 

The  simulation  assumes  the  test  vehicle  possesses  the  same 
configuration  of  RCS  jets  as  the  actual  AFE  vehicle.  The  simulation 
uses  the  jet  properties  listed  in  Table  6.3. 

7.2.3  Equations  of  Motion 

As  the  aerocapture  vehicle  descends  deeper  into  the 
atmosphere,  the  aerodynamic  forces  and  torques  experienced  by  the 
vehicle  increase.  The  aerodynamics  are  simulated  using  the 
aerodynamic  properties  of  the  actual  AFE  vehicle  [37].  A  simple 
exponential  atmospheric  model  based  on  the  1962  U.S.  Standard 
Atmosphere  [6]  is  used  to  provide  the  atmospheric  density 
encountered  for  each  given  altitude. 

The  overall  simulation  is  cycled  at  25  hz.  The  translational  and 
rotational  states  of  the  aerocapture  vehicle  are  specified  at  each  step 
by  the  12-element  state  vector 


Xz  -  distance  from  the  center  of  Earth 
|i  =  latitude 
^  =  longitude 

VI  =  North  component  of  velocity 
V2  =  East  component  of  velocity 

V3  =  radial  (toward  center  of  Earth)  velocity  component 
6,  (p,  Y  =  Euler  angles  relating  the  body  axes  to  the 

velocity  vector  as  described  in  Appendix  A 

p  =  roll  rate 
q  =  pitch  rate 
r  =  yaw  rate 


During  each  step,  a  first  order  differential  equation  for  each  of  the 
elements  in  the  state  vector  is  calculated  based  on  the  jets  fired  and 
the  aerodynamic  forces  and  torques  encountered.  These  equations  of 
motion  are  then  integrated  using  a  fourth-order  Runge-Kutta  algo¬ 
rithm  [38]  to  determine  the  new  states  at  the  end  of  each  step.  The 
current  velocity  angles  are  then  calculated  based  on  the  velocity 
elements  (i.e.,  vi,  V2,  and  V3)  and  the  Euler  angles. 

7.3  Robustness  Testing 

The  main  objective  of  the  candidate  design  was  to  integrate  a 
mass  property  identification  algorithm  with  a  robust  control  law  and 
an  adaptable  jet  selection  algorithm  to  extend  the  range  of  mass 
properties  for  which  the  attitude  control  system  is  robust.  To 
demonstrate  the  improvement  of  the  new  design,  robustness  testing 
was  conducted  for  the  following  three  control  system  configurations: 

1.  AFE  DAP  [16] 
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2. 


new  control  law/jet  selection  algorithm  with  mass 
property  estimates  set  to  the  nominal  values  (i.e,,  without 
the  mass  property  identification  algorithm) 

3.  complete  candidate  attitude  control  system  design  (i.e., 
employing  the  mass  property  identification  algorithm) 

Tests  with  the  AFE  DAP  are  performed  only  for  comparison 
purposes.  The  range  of  mass  properties  to  which  the  AFE  DAP  is 
subjected  during  these  tests  will  not  be  encountered  during  the 
actual  AFE  mission.  The  AFE  DAP  is  not  designed  nor  required  to  be 
capable  of  providing  attitude  control  over  this  range  of  mass 
properties. 

For  each  c  ..figuration,  a  60  sec  simulation  is  repeatedly 

performed  with  the  mass  properties  changed  from  one  simulation 

run  to  the  next.  Figure  7.1  shows  the  results  of  a  single  run  of  the  60 

sec  simulation  for  the  complete  candidate  control  system  (i.e., 

configuration  #3)  and  a  vehicle  with  the  nominal  mass  properties. 
The  commanded  bank  angle  is  initially  45  deg  and  then  varies 
sinusoidally.  The  aerocapture  vehicle  begins  with  the  initial 
conditions  listed  in  Section  7.2.1,  reaches  the  commanded  trajectory 
in  approximately  8  sec,  and  then  begins  tracking  while  maintaining 
the  commanded  a  and  p. 

Nominal  mass  properties  are  used  during  the  first  run  of  the  60 
sec  simulation  for  each  configuration.  For  each  successive  run,  the 
center  of  mass  location  is  moved  1  inch  from  the  previous  location 
(e.g.,  on  the  10*^  run,  the  center  of  mass  is  moved  9  inches  from  the 
nominal  position).  Each  1  inch  move  is  achieved  by  moving  the 
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center  of  mass  an  equal  distance  along  the  +x,  +y,  and  +z  axes  (i.e., 
0.58  inches  along  each  axis).  The  inertia  matrix  is  changed  as 
described  in  Appendix  C.  This  is  done  to  simulate  changes  in  both 
the  center  of  mass  location  and  the  inertia  matrix,  the  mass 
properties  that  affect  performance  of  the  control  system. 


Figure  7.1.  Velocity  Angle  Response  of  Typical  60  Sec 

Simulation 


Figures  7.2  and  7.3  show  the  testing  results  for  the  three 
control  system  configurations.  Figure  7.2  is  a  plot  of  the  total 
number  of  jet  firings  required  during  each  60  sec  simulation  vs.  the 
displacement  of  the  center  of  mass  location  for  each  simulation. 
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Figure  7.3  shows  the  maximum  bank  angle  tracking  error,  after  the 
desired  trajectory  is  reached,  for  each  simulation. 


center  of  mass  displacement  (inches) 

Figure  7.2.  Jet  Firings  vs.  Center  of  Mass  Displacement 

For  the  AFE  DAP,  the  number  of  jet  firings  begins  to  increase 
rapidly  when  the  center  of  mass  location  is  displaced  by  more  than 
33  inches  from  the  nominal  location.  The  bank  angle  errors  increase 
dramatically  after  about  40  inches  of  displacement.  For  the  new 
design  without  the  mass  property  identification  algorithm  (i.e., 
configuration  #2),  the  number  of  jet  firings  increases  rapidly  at  44 
inches  of  displacement  while  the  maximum  bank  angle  errors  begin 
to  increase  rapidly  after  about  60  inches  of  displacement.  For  the 
new  design  employing  the  mass  property  identification  algorithm 
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(i.e.,  configuration  #3),  bank  angle  errors  do  not  increase  significantly 
over  the  range  of  mass  properties  tested.  The  number  of  jet  firings 
does  slowly  increase,  however,  as  the  center  of  mass  is  displaced 
further  from  the  nominal  location.  This  is  expected  since  the 
magnitudes  of  the  inertia  matrix  elements  increase  as  the  center  of 
mass  location  is  displaced  further,  requiring  additional  control 
activity  to  move  the  larger  inertia. 


center  of  mass  displacement  (inches) 

Figure  7.3.  Max  Bank  Angle  Error  vs.  Center  of  Mass 

Displacement 


The  displacements  at  which  the  bank  angle  errors  rapidly 
increase  can  be  considered  the  points  at  which  the  control  system 
becomes  unstable.  The  AFE  DAP  becomes  unstable  at  a  40  inch 
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displacement  of  the  center  of  mass  location.  The  new  design  using 
the  nominal  values  for  mass  property  estimates  becomes  unstable  at 
60  inches  of  displacement,  more  robust  than  the  AFE  DAP,  but  not 
robust  to  the  full  range  of  mass  properties  tested.  The  complete 
design  remains  stable  over  the  full  range  of  mass  properties, 
demonstrating  the  benefit  of  incorporating  the  mass  property 
identification  algorithm. 

7.4  Performanc«»  of  the  Candidate  Attitude  Control 
System 

Section  7.3  demonstrated  the  benefit  of  using  a  mass  property 
identification  algorithm  to  extend  the  range  of  mass  properties  for 
which  a  control  system  is  robust.  For  each  run  of  the  60  sec 
simulation,  the  actual  mass  properties  were  set  as  described  above 
but  remained  constant  during  the  simulation.  The  identification 
algorithm  was  able  to  accurately  estimate  the  mass  properties  at  the 
beginning  of  each  simulation,  and  the  control  system  was  able  to 
track  the  simple  desired  trajectory. 

This  section  demonstrates  the  ability  of  the  candidate  control 
system  to  track  a  bank  angle  trajectory  that  includes  a  simulated 
commanded  bank  reversal,  along  with  the  ability  to  detect,  identify, 
and  adapt  to  mass  properties  that  change  during  the  simulation. 
Aerodynamic  torques  are  not  simulated  during  the  tests  described  in 
this  section,  but  are  simulated  during  the  tests  presented  in  Section 
7.5.  Section  7.4.1  describes  the  performance  of  the  candidate  control 
system  during  tracking  of  the  desired  bank  angle  trajectory  when 
the  vehicle  has  the  nominal  mass  properties.  In  Section  7.4.2,  the 
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candidate  control  system  tracks  the  desired  trajectory  while  also 
detecting  and  adapting  to  large  step  changes  in  the  vehicle  mass 
properties.  Finally,  Section  7.4.3  shows  the  ability  of  the  candidate 
control  system  to  detect  and  adapt  to  gradual  mass  property  changes. 

The  remaining  tests  use  the  same  commanded  bank  angle 
trajectory.  This  trajectory  was  designed  to  mimic  the  main 
characteristics  of  an  aerocapture  bank  angle  trajectory  commanded 
by  the  AFE  guidance  law.  To  provide  maximum  control  margin,  the 
in-plane  portion  of  the  AFE  guidance  commands  most  of  the  lift  out 
of  the  desired  vertical  plane.  If  the  aerocapture  vehicle  slows  too 
quickly,  guidance  then  can  command  more  in-plane  lift  to  pull  the 
vehicle  out  of  the  atmosphere.  If  the  vehicle  does  not  decelerate 
quickly  enough,  more  in-plane  lift  can  be  commanded  to  pull  the 
vehicle  deeper  into  the  atmosphere  increasing  deceleration.  Once  the 
out-of-plane  lift  pulls  the  aerocapture  vehicle  out  of  the  desired 
vertical  plane,  a  bank  reversal  is  commanded  by  the  out-of-plane 
guidance  to  bring  the  vehicle  back  to  the  desired  plane. 

To  simulate  these  characteristics,  the  desired  bank  angle 
trajectory  calls  for  an  initial  bank  angle  of  75  deg.  It  then  varies  30 
deg  sinusoidally  around  90  deg  bank  until  60  sec  into  the  simulation. 
At  60  sec,  the  mid  point  of  the  simulation,  a  bank  angle  reversal  is 
commanded.  The  final  60  sec  of  the  simulation  calls  for  a 
sinusoidally  varying  trajectory  centered  at  -90  deg  bank. 

7.4.1  Nominal  Performance 

Figure  7.4  illustrates  the  tracking  performance  of  the  candidate 
aerocapture  attitude  control  system  for  a  vehicle  with  the  initial 
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conditions  and  the  nominal,  unchanging  mass  properties  listed  in 
Section  7.2.  Even  with  the  initial  0.5  sec  mass  property  identification 
period  when  the  control  law  and  jet  selection  algorithm  are 
bypassed,  the  desired  bank  angle  trajectory  is  reached  in  under  13 
sec,  in  agreement  with  Equation  3.12 

13.0  sec 

Tli  0.1 


time  (sec) 

Figure  7.4.  Velocity  Angle  Response  -  New  Design 
(Nominal  Mass  Properties) 


After  the  initial  reach  period,  bank  angle  tracking  precision  is 
within  the  desired  2  deg  until  the  commanded  bank  reversal  when  a 
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new  reaching  period  is  begun.  During  the  entire  simulation,  a  and  P 
are  precisely  maintained  at  their  commanded  values. 

Figure  7.5  shows  the  value  of  the  sliding  variable  si  associated 
with  the  bank  angle  during  the  simulation.  Initially,  and  at  the  point 
of  the  commanded  bank  reversal,  the  value  of  si  is  off  the  sliding 
surface  defined  by  si(x,t)=0.  The  sliding  surface  is  reached  and 
maintained  to  within  the  boundary  layer  0 1  in  the  time  specified  by 
Equation  3.12. 


Figure  7.5.  Sliding  Variable  and  Boundary  Layer 
(Nominal  Mass  Properties) 
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Table  7,1  shows  additional  results  displayed  to  the  screen 
during  the  simulation.  At  t=0.52  sec  into  the  simulation,  the 
estimated  mass  properties  are  reported  and  supplied  to  the  control 
law  and  the  jet  selection  algorithm.  The  mass  property  estimates  are 
accurate  to  within  1%  of  the  actual  values.  Table  7.1  also  shows  how 
often  each  jet  is  fired  during  the  simulation.  For  the  entire 
simulation,  554  jet  firings,  each  of  40  msec  duration,  are  performed 
to  provide  attitude  control,  while  12  jet  firings,  also  of  40  msec 
duration,  are  used  for  mass  property  estimation. 


Table  7.1  New  Design  (Nominal  Simulation) 

***  TIME  =  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2243.93  -4.97  -63.49 

-4.97  1459.20  21.03 

-63.49  21.03  1178.62 

ESTIMATED  CENTER  OF  MASS  IS: 

83.66 

-.03 

-.61 
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Table 

7.1  (continued) 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  CONTROL  LAW/ JET  SELECTION; 

JET  #  1 

99  JET  #  9 

0 

JET  #  2 

137  JET  #10 

0 

JET  #  3 

111  JET  #11 

0 

JET  #  4 

114  JET  #12 

0 

JET  #  5 

59  JET  #13 

0 

JET  #  6 

1  JET  #14 

0 

JET  #  7 

10  JET  #15 

0 

JET  #  8 

23  JET  #16 

0 

TOTAL  #  JET  FIRINGS (CONTROL 

LAW/ JET  SELECT) :  554 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID: 

JET  #  1 

1  JET  #  9 

0 

JET  #  2 

0  JET  #10 

0 

JET  #  3 

0  JET  #11 

0 

JET  #  4 

0  JET  #12 

0 

JET  #  5 

4  JET  #13 

0 

JF'f  #  6 

4  JET  #14 

0 

JET  #  7 

2  JET  #15 

0 

JET  #  8 

1  JET  #16 

0 

TOTAL  #  JET  FIRINGS (MASS  PROP  ID) :  12 

For  comparison  purposes,  the  same  simulation  is  performed 
with  the  AFE  DAP.  Figure  7.6  shows  the  tracking  response  to  be 
similar  to  that  for  the  candidate  control  system,  but  somewhat  less 
precise.  Table  7.2  indicates  that  the  AFE  DAP  uses  a  total  of  825  jet 
firings  for  attitude  control  for  this  simulation,  nearly  50%  more  jet 
firings  than  the  new  design  requires. 


142 


velocity  angles  (deg) 


Figure  7.6.  Velocity  Angle  Response  -  AFE  DAP 
(Nominal  Mass  Properties) 
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Table  7.2  AFE  DAP  (Nominal  Simulation) 

***  TIME  -  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 


NUMBER  OF  JET  FIRINGS  COMMANDED  BY  CONTROL  LAW/ JET  SELECTION: 


JET 

# 

1 

165 

JET 

#  9 

JET 

« 

2 

165 

JET 

#10 

JET 

# 

3 

159 

JET 

#11 

JET 

* 

4 

159 

JET 

#12 

JKT 

# 

5 

81 

JET 

#13 

JET 

# 

6 

2 

JET 

#14 

JET 

# 

7 

52 

JET 

#15 

JET 

# 

8 

42 

JET 

#16 

TOTAL  #  JET  FIRINGS (CONTROL  LAW/JET  SELECT):  825 


7.4.2  Performance  in  the  Presence  of  Large  Step  Changes  to 
the  Mass  Properties 

The  above  simulation  was  repeated  with  step  changes  in  the 
mass  properties  to  demonstrate  the  ability  of  the  candidate  control 
system  to  detect  and  adapt  to  large  changes  in  the  vehicle  mass 
properties.  The  simulation  begins  with  the  vehicle  having  the 
nominal  mass  properties  listed  in  Table  6.2.  At  t=30  sec,  the  center 
of  mass  is  displaced  45  inches  from  its  nominal  location  with  the  new 
inertia  matrix  calculated  as  in  Appendix  C.  At  t=90  sec,  another  step 
change  is  simulated  by  changing  the  mass  properties  back  to  the 
original  nominal  values. 
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Figure  7.7.  Velocity  Angle  Response  -  New  Design 
(Large  Step  Changes  in  Mass  Properties) 

Figure  7.7  shows  the  tracking  response  of  the  candidate  control 
system.  Comparing  to  Figure  7.4,  little  change  in  tracking 
performance  is  discernable.  Table  7.3  lists  the  significant  results  of 
this  test.  The  initial  mass  property  estimates  at  t=0.52  sec  are  the 
same  as  the  initial  estimates  found  in  Table  7.1  for  the  nominal 
simulation.  At  t=32.08  sec,  the  first  step  change  is  detected.  After 
14  jet  firings  commanded  by  the  special  input  selection  scheme  in 
the  mass  property  identification  algorithm,  the  mass  property 
estimates  are  updated  and  supplied  to  the  control  law  and  jet 
selection  algorithm  at  t=34.92  sec.  The  second  mass  property  change 
is  detected  at  t=91.72  sec.  After  5  jet  firings  commanded  by  the 
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mass  property  identification  algorithm,  the  identification  period 
ends,  and  an  accurate  update  to  the  mass  property  estimates  occurs 
at  t=92.76  sec. 


Table  7.3  New  Design  (Step  Changes  in  Mass 

Properties) 

***  TIME  =  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2243,93  -4.97  -63.49 

-4.97  1459.20  21.03 

-63.49  21.03  1178.62 

ESTIMATED  CENTER  OF  MASS  IS: 

83.66 

-.03 

-.61 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  32.08  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 


JET 

# 

7 

AT 

T 

= 

32.20 

SEC 

JET 

« 

5 

AT 

T 

32.40 

SEC 

JET 

« 

5 

AT 

T 

= 

32.60 

SEC 

JET 

# 

7 

AT 

T 

= 

32.80 

SEC 

JET 

* 

2 

AT 

T 

= 

33.00 

SEC 

JET 

# 

6 

AT 

T 

33.20 

SEC 

JET 

# 

5 

AT 

T 

= 

33.40 

SEC 

JET 

« 

1 

AT 

T 

= 

33.60 

SEC 

JET 

« 

5 

AT 

T 

33.80 

SEC 

JET 

* 

4 

AT 

T 

= 

34.00 

SEC 

JET 

» 

2 

AT 

T 

34.20 

SEC 

JET 

# 

5 

AT 

T 

= 

34 . 40 

SEC 

JET 

# 

6 

AT 

T 

= 

34.60 

SEC 

JET 

# 

5 

AT 

T 

= 

34.80 

SEC 

Table  7.3  (continued) 


***  TIME  =  34.92  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

4554.31  -1158.06  -1216.66 
-1158.06  3765.31  -1132.06 

-1216.66  -1132.06  3484.71 

CENTER  OF  MASS  IS: 

109.58 

25.95 

25.37 

UPDATED  MASS  PROPERTY  ESTIMATES: 
EST  INERTIA  MATRIX  IS: 

4554.30  -1158.20  -1216.74 
-1158.20  3765.62  -1132.00 

-1216.74  -1132.00  3484.70 

ESTIMATED  CENTER  OF  MASS  IS: 
109.58 
25.95 
25.37 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  91.72  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 
JET  #  5  AT  T  =  91.80  SEC 

JET  #  7  AT  T  =  92.00  SEC 

JET  #  8  AT  T  =  92.20  SEC 

JET  #  8  AT  T  =  92.40  SEC 

JET  #  5  AT  T  =  92.60  SEC 

***  TIME  =  92.76  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2248.21  -4.97  -63.58 

-4.97  1459.02  21.05 

-63.58  21.05  1178.57 

ESTIMATED  CENTER  OF  MASS  13: 

83.60 

-.03 

-.61 
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Table 

7.3  (continued) 

NXniBER  OF  JET  FIRINGS  COMMANDED  BY  CONTROL  LAW/ JET  SELECTION: 

JET  #  1 

315  JET  #  9 

0 

JET  #  2 

108  JET  #10 

0 

JET  #  3 

287  JET  #11 

0 

JET  #  4 

85  JET  #12 

0 

JET  #  5 

55  JET  #13 

0 

JET  #  6 

27  JET  #14 

0 

JET  #  7 

56  JET  #15 

0 

JET  #  8 

63  JET  #16 

1 

1 

1  TOTAL  #  JET  FIRINGS (CONTROL 

LAW/JET  SELECT) :  997 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID: 

JET  #  1 

2  JET  #  9 

0 

JET  #  2 

2  JET  #10 

0 

JET  #  3 

0  JET  #11 

0 

JET  #  4 

1  JET  #12 

0 

JET  #  5 

11  JET  #13 

0 

JET  #  6 

6  JET  #14 

0 

JET  #  7 

5  JET  #15 

0 

JET  #  8 

3  JET  #16 

0 

TOTAL  #  JET  FIRINGS (MASS  PROP  ID) :  30 

For  comparison,  the  AFE  DAP  is  run  against  this  simulation  and 
the  tracking  results  are  shown  in  Figure  7.8.  The  initial  tracking 
response  is  good,  but  the  vehicle  becomes  unstable  when  the  AFE 
DAP  attempts  to  implement  the  bank  reversal  after  the  first  step 
change  in  the  mass  properties.  Even  after  the  mass  properties  return 
to  iheir  nominal  values  at  t=90  sec  the  AFE  DAP  cannot  recover 
before  the  end  of  the  simulation. 
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Figure  7.8.  Velocity  Angle  Response  -  AFE  DAP 
(Large  Step  Changes  in  Mass  Properties) 


7.4.3  Performance  in  the  Presence  of  Gradually  Changing 
Mass  Properties 

To  demonstrate  the  ability  of  the  candidate  control  system  to 
identify  and  adapt  to  gradual  changes  in  mass  properties,  the 
simulation  was  again  repeated  but  with  the  center  of  mass 
continually  moving.  Beginning  with  the  nominal  mass  properties,  the 
center  of  mass  is  moved  0.5  inches/sec  throughout  the  full 
simulation.  As  in  the  testing  above,  each  0.5  inch  move  is  achieved 
by  moving  the  center  of  mass  an  equal  distance  along  the  +x,  +y,  and 
+z  axes.  At  the  end  of  the  120  sec  simulation,  the  center  of  mass 
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velocity  angles  (deg) 


location  is  displaced  a  total  of  60  inches  from  the  nominal  location  at 
the  coordinates 


rcm(l20  sec)  = 


118.24 
34.61 
.  34.03  . 


inches 


3) 


Throughout  the  simulation,  the  inertia  matrix  remains  equal  to  the 
nominal  inertia  matrix. 


Figure  7.9.  Velocity  Angle  Response  -  New  Design 
(Slow,  Continual  Mass  Property  Changes) 
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Figure  7.9  shows  the  tracking  response  of  the  candidate  control 
system  for  this  test.  Again  comparing  to  Figures  7.4  and  7.7,  no 
change  in  tracking  performance  is  discernable.  Table  7.4  lists  the 
main  events  during  this  test.  The  initial  mass  property  estimates  are 
supplied  by  the  identification  algorithm  at  t=0.52  sec.  At  t=37.32  sec, 
the  first  detection  of  a  mass  property  change  occurs.  Mass  property 
identification  is  performed  while  control  is  maintained,  and  accurate, 
updated  mass  property  estimates  are  supplied  at  t=38.68  sec.  The 
center  of  mass  continues  to  be  moved,  and  at  t=69.92  sec,  the  flag  is 
set  indicating  another  mass  property  change  has  been  detected.  The 
new  identification  period  lasts  until  t=71.16  sec  when  the  new 
estimates  are  supplied.  A  final  change  is  detected  at  t=109.52  sec 
with  updated  mass  property  estimates  supplied  at  t=l  10.80  sec.  The 
simulation  ends  before  another  mass  property  change  is  detected. 


Table 

7.4  New  Design 

(Slow,  Continual 

Changes  in  Mass  Properties) 

***  TIME  =  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00 

-63.60 

-5.00  1459.20 

21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.75 

.12 

-.46 

Table  7.4  (continued) 


UPDATED  MASS  PROPERTY  ESTIMATES: 
EST  INERTIA  MATRIX  IS: 

2243.42  -4.97  -63.48 

-4.97  1459.14  21.07 

-63.48  21.07  1178.54 

ESTIMATED  CENTER  OF  MASS  IS: 
83.66 
.04 
-.54 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  37.32  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 
JET  #  7  AT  T  =  37.40  SEC 

JET  #  1  AT  T  =  37.60  SEC 

JET  #  2  AT  T  =  37.80  SEC 

JET  #  5  AT  T  =  38.00  SEC 

JET  #  7  AT  T  =  38.20  SEC 

JET  #  3  AT  T  =  38.40  SEC 

JET  #  6  AT  T  =  38.60  SEC 

***  TIME  =  38.68  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

94,77 

11.14 

10.56 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2247.74  -4.78  -63.59 

-4.78  1457.49  21,56 

-63.59  21.56  1178.72 

ESTIMATED  CENTER  OF  MASS  IS: 

94.62 

11.04 

10.46 
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Table  7.4  (continued) 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  69.92  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 
JET  #  7  AT  T  =  70.00  SEC 

JET  #  7  AT  T  =  70.20  SEC 

JET  #  6  AT  T  =  70.40  SEC 

JET  #  5  AT  T  =  70.60  SEC 

JET  #  6  AT  T  =  70.80  SEC 

JET  #  4  AT  T  =  71.00  SEC 

***  TIME  -  71.16  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

104.14 

20.51 

19.93 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2265.55  2.38  -58.89 

2.38  1395.09  30.76 

-58.89  30.76  1119.01 

ESTIMATED  CENTER  OF  MASS  IS: 

101.80 

20.29 

20.51 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  109.52  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 
JET  #  7  AT  T  =  109.60  SEC 

JET  #  1  AT  T  =  109.80  SEC 

JET  #  2  AT  T  =  110.00  SEC 

JET  #  2  AT  T  =  110.20  SEC 

JET  #  6  AT  T  =  110.40  SEC 

JET  #  5  AT  T  -  110.60  SEC 


Table  7.4  (continued) 


***  TIME  -  110.80  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

115.59 

31.96 

31.38 

UPDATED  MASS  PROPERTY  ESTIMATES: 
EST  INERTIA  MATRIX  IS: 

2247.80  -4.05  -63.83 

-4.05  1458.48  20.70 

-63.83  20.70  1176.29 

ESTIMATED  CENTER  OF  MASS  IS: 

115.46 

31.79 

31.29 


NUMBER  OF  JET  FIRINGS  COMMANDED  BY  CONTROL  LAW/ JET  SELECTION: 


JET  #  1 

131  JET  #  9 

0 

JET  #  2 

125  JET  #10 

0 

JET  #  3 

133  JET  #11 

0 

JET  #  4 

108  JET  #12 

0 

JET  #  5 

87  JET  #13 

0 

JET  #  6 

2  JET  #14 

0 

JET  #  7 

8  JET  #15 

0 

JET  #  8 

51  JET  #16 

0 

TOTAL  #  JET  FIRINGS (CONTROL  LAW/ JET  SELECT):  645 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID: 


JET  #  1 

3  JET  #  9 

0 

JET  #  2 

3  JET  #10 

0 

JET  #  3 

1  JET  #11 

0 

JET  #  4 

1  JET  #12 

0 

JET  #  5 

7  JET  #13 

0 

JET  #  6 

8  JET  #14 

0 

JET  #  7 

7  JET  #15 

0 

JET  #  8 

1  JET  #16 

0 

TOTAL  #  JET  FIRINGS (MASS  PROP  ID) :  31 
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time  (sec) 

Figure  7.10.  Velocity  Angle  Response  -  AFE  DAP 
(Slow,  Continual  Mass  Property  Changes) 

Figure  7.10  shows  the  tracking  response  of  the  AFE  DAP  run 
against  this  simulation.  At  the  midpoint  (i.e.,  t=60  sec)  the  center  of 
mass  has  only  been  moved  30  inches.  The  AFE  DAP  is  still  stable 
with  this  amount  of  displacement  and  is  able  to  implement  the  bank 
reversal.  At  roughly  t=100  sec,  the  limits  of  stability  have  been 
exceeded,  and  the  system  becomes  unstable. 

7.5  Effects  of  Aerodynamics 

In  Section  1.5.2  the  assumption  was  made  that  the 
aerodynamic  torques  experienced  by  the  vehicle  during  the 
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aerocapture  maneuver  are  small  compared  to  the  torques  applied  by 
the  RCS  jets.  The  control  law  was  designed  to  be  robust  to  the 

aerodynamic  torques,  and  the  candidate  attitude  control  system 
should  be  robust  to  the  disturbances  caused  by  the  aerodynamic 
torques.  The  tests  performed  in  Section  7.4  were  repeated  with  the 
aerodynamic  torques  simulated  to  verify  these  assumptions.  The 

results  presented  below  show  that  the  control  system  is  indeed 
robust  to  these  aerodynamic  disturbances.  The  results  also  show, 
however,  that  the  mass  property  identification  algorithm  has 
difficulty  in  converging  to  accurate  estimates  when  the  vehicle  is  at 
its  deepest  points  in  the  atmosphere  where  the  magnitudes  of  the 
aerodynamic  torques  are  at  their  greatest.  This  is  not  surprising 

since  in  Section  5.2  for  the  filter  design  for  mass  property  estimation, 
the  additional  assumption  was  made  that  the  aerodynamic  torques 
were  small  enough  to  be  treated  as  noise  to  the  rate  change 

measurements.  The  mass  property  identification  algorithm  design 
has  not  yet  been  extended  to  accurately  estimate  mass  properties  in 
the  presence  of  significant  aerodynamic  disturbances. 

7.5.1  Nominal  Performance 

For  a  vehicle  with  the  nominal,  unchanging  mass  properties, 
the  response  of  the  candidate  control  system  is  very  similar  to  the 
response  described  in  Section  7.4.1.  The  tracking  response  is  almost 
identical,  and  Table  7.5  shows  that  the  initial  mass  property 
estimates  are  nearly  the  same  as  the  initial  estimates  listed  in  Table 
7.1.  The  number  of  jet  firings  commanded  by  the  jet  selection  has 
risen  slightly,  from  554  to  586,  to  counter  the  aerodynamic  torques. 
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Table  7.5  New  Design  w/  Aerodynamics 
(Nominal  Simulation) 

***  TIME  =  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2243.57  -4.97  -63.47 

-4.97  1459.19  21.04 

-63.47  21.04  1178.61 

ESTIMATED  CENTER  OF  MASS  IS: 

83.66 

-.03 

-.62 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  CONTROL  LAW/ JET  SELECTION: 


JET  #  1 

97  JET  #  9 

0 

JET  #  2 

145  JET  #10 

0 

JET  #  3 

120  JET  #11 

0 

JET  #  4 

111  JET  #12 

0 

JET  #  5 

84  JET  #13 

0 

JET  #  6 

1  JET  #14 

0 

JET  #  7 

6  JET  #15 

0 

JET  #  8 

22  JET  #16 

0 

TOTAL  #  JET  FIRINGS (CONTROL  LAW/JET  SELECT):  586 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID: 


JET  #  1 

1  JET  #  9 

0 

JET  #  2 

0  JET  #10 

0 

JET  #  3 

0  JET  #11 

0 

JET  #  4 

0  JET  #12 

0 

JET  #  5 

4  JET  #13 

0 

JET  #  6 

4  JET  #14 

0 

JET  #  7 

2  JET  #15 

0 

JET  #  8 

1  JET  #16 

0 

TOTAL  #  JET  FIRINGS (MASS  PROP  ID):  12 


7.5.2  Performance  in  the  Presence  of  Large  Step  Changes  to 
the  Mass  Properties 

Figure  7.11  shows  the  altitude  of  the  vehicle  throughout  each 
simulation. 


Figure  7.11.  Altitude  vs.  Time 

The  aerodynamic  torques  do  have  a  significant  impact  during 
the  test  in  which  large  step  changes  in  the  mass  properties  are 
simulated.  Table  7.6  shows  that  the  candidate  system  is  still  able  to 
identify  and  adapt  to  the  first  change  at  t=30  sec.  The  change  is 
detected  at  t=33.76  sec  and  the  updated  mass  property  estimates  are 
supplied  at  t=36.08  sec.  These  estimates,  however,  are  not  as 
accurate  as  the  first  updated  estimates  listed  in  Table  7.3.  The 
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aerodynamic  torques  become  larger  as  the  vehicle  descends  through 
the  atmosphere  and  begin  degrading  the  accuracy  of  the  mass 
property  estimates. 


Table  7.6  New  Design  w/  Aerodynamics 
(Step  Changes  in  Mass  Properties) 

***  TIME  -  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2243,57  -4,97  -63.47 

-4.97  1459.19  21,04 

-63.47  21.04  1178.61 

ESTIMATED  CENTER  OF  MASS  IS: 

83.66 

-.03. 

-.62 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  33.76  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 


JET 

# 

7 

AT 

T 

= 

33.80 

SEC 

JET 

# 

5 

AT 

T 

= 

34.00 

SEC 

JET 

# 

3 

AT 

T 

= 

34.20 

SEC 

JET 

# 

7 

AT 

T 

= 

34.40 

SEC 

JET 

# 

8 

AT 

T 

= 

34.60 

SEC 

JET 

# 

4 

AT 

T 

= 

34.80 

SEC 

JET 

# 

5 

AT 

T 

= 

35.00 

SEC 

JET 

» 

3 

AT 

T 

- 

35.20 

SEC 

JET 

« 

8 

AT 

T 

35.40 

SEC 

JET 

# 

4 

AT 

T 

= 

35.60 

SEC 

JET 

« 

7 

AT 

T 

- 

35.80 

SEC 

JET 

# 

2 

AT 

T 

= 

36.00 

SEC 
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Table 

7.6 

New 

Design  w/  Aerodynamics 

(Step  Changes  in  Mass  Properties) 

TIME  =  36.08  SEC 

*  *  * 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX 

IS: 

4554,31  - 

1158, 

06  -1216.66 

-1158.06 

3765. 

31  -1132.06 

-1216.66  - 

1132. 

06  3484.71 

CENTER  OF  MASS 

IS: 

109.58 

25.95 

25.37 

UPDATED  MASS  PROPERT^^ 

ESTIMATES: 

EST  INERTIA  MATRIX 

IS: 

4506.64  - 

1135. 

88  -1205.18 

-1135,88 

3779. 

37  -1138.18 

-1205.18  - 

1138. 

18  3467.13 

ESTIMATED  CENTER  OF 

MASS 

IS: 

109.49 

24.55 

24.97 

MASS  PROPERTY  CHANGE 

DETECTED  AT  T  «  90.36  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 

JET  #  5  AT 

T  - 

90.40 

SEC 

JET  #  7  AT 

T  = 

90.60 

SEC 

JET  #  5  AT 

T  = 

90.80 

SEC 

.  JET  #  5  AT 

T  = 

91.00 

SEC 

JET  #  5  AT 

T  = 

91.20 

SEC 

JET  #  6  AT 

T  = 

91.40 

SEC 

JET  #  8  AT 

T  = 

91.60 

SEC 

JET  #  5  AT 

T  - 

91.80 

SEC 

JET  #  5  AT 

T  = 

92.00 

SEC 

JET  #  6  AT 

T  - 

92.20 

SEC 

JET  #  8  AT 

T  = 

92.40 

SEC 

JET  #  8  AT 

T  = 

92.60 

SEC 

JET  #  8  AT 

T  - 

92.80 

SEC 

JET  #  8  AT 

T  = 

S3  00 

SEC 

JET  #  8  AT 

T  = 

93.20 

SEC 

JET  #  5  AT 

T  = 

93.40 

SEC 

JET  #  5  AT 

T  - 

93.60 

SEC 

JET  #  5  AT 

T  - 

93.80 

SEC 

JET  #  5  AT 

T  - 

94.00 

SEC 

JET  #  8  AT 

T  = 

94.20 

SEC 

JET  #  8  AT 

T  - 

94.40 

SEC 

JET  #  5  AT 

T  ” 

94.60 

SEC 

JET  #  8  AT 

T  = 

94.80 

SEC 

JET  #  5  AT 

T  = 

95.00 

SEC 
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Table  7.6  (continued) 


JET 

# 

5 

AT 

T 

= 

95.20 

SEC 

JET 

# 

8 

AT 

T 

95.40 

SEC 

JET 

# 

5 

AT 

T 

= 

95.60 

SEC 

JET 

# 

5 

AT 

T 

= 

95.80 

SEC 

JET 

# 

8 

AT 

T 

= 

96.00 

SEC 

JET 

# 

5 

AT 

T 

= 

96.20 

SEC 

JET 

# 

8 

AT 

T 

= 

96.40 

SEC 

JET 

# 

8 

AT 

T 

- 

96.60 

SEC 

JET 

# 

5 

AT 

T 

= 

96.80 

SEC 

JET 

# 

8 

AT 

T 

= 

97.00 

SEC 

JET 

# 

5 

AT 

T 

= 

97.20 

SEC 

JET 

# 

4 

AT 

T 

= 

97.40 

SEC 

JET 

# 

5 

AT 

T 

- 

97.60 

SEC 

JET 

# 

8 

AT 

T 

= 

97.80 

SEC 

JET 

# 

5 

AT 

T 

- 

98.00 

SEC 

***  TIME  -  98.20  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.60 

-.03 

-.61 

UPDATED  MASS  PROPERTY  ESTIMATES: 
EST  INERTIA  MATRIX  IS: 

2481.96  34.35  -1.22 

34.35  1414.26  1.85 

-1.22  1.85  1146.97 

ESTIMATED  CENTER  OF  MASS  IS: 
79.67 
3.46 
3.22 


NUMBER  OF  JET  FIRINGS  COMMANDED 
JET  #  1:  360  JET  #  9: 

JET  #  2:  110  JET  #10: 

JET  #3:  333  JET  #11: 

JET  #4:  75  JET  #12: 

JET  #5:  85  JET  #13: 

JET  #6:  47  JET  #14: 

JET  #7:  57  JET  #15: 

JET  #8:  61  JET  #16: 


BY  CONTROL  LAW/JET  SELECTION 
0 
0 
0 
0 
0 
0 
0 
0 


TOTAL  #  JET  FIRINGS (CONTROL  LAW/JET  SELECT):  1128 


Table 

7.6  (continued) 

NUMBER  OF  JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID: 

JET  #  1 

1  JET  #  9 

0 

JET  #  2 

1  JET  #10 

0 

JET  #  3 

2  JET  #11 

0 

JET  #  4 

3  JET  #12 

0 

JET  #  5 

23  JET  #13 

0 

JET  #  6 

5  JET  #14 

0 

JET  #  7 

6  JET  #15 

0 

JET  #  8 

18  JET  #16 

0 

TOTAL  #  JET  FIRINGS (MASS  PROP  ID) :  59 

Table  7.6  also  shows  that  the  second  mass  property  change  is 
detected  at  t=90.36  sec.  This  time  numerous  jet  firings  are  required 
before  the  filter  can  converge  to  mass  property  estimates  accurate 
enough  to  end  the  identification  period.  After  the  first  few  jet 
firings,  the  input  selection  scheme  in  the  mass  property  identification 
algorithm  commands  mainly  jets  #5  and  #8.  This  is  an  indication 
that  the  filter  has  become  optimistic  (i.e.,  the  trace  elements  of  P 
have  decreased  nearly  to  their  values  before  the  identification  period 
began).  Little  information  is  extracted  from  these  later  jet  firings 

and  the  filter  has  difficulty  in  converging  to  accurate  estimates. 

Finally,  at  t=98.20  sec,  the  two  criteria  described  in  Section  6.4.3  are 
met  and  updated  mass  property  estimates  are  supplied.  Though  the 
accuracy  has  been  degraded  by  the  aerodynamic  torques,  the 

estimates  are  accurate  enough  for  the  control  system  to  track  the 

desired  trajectory. 
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7.5.3  Performance  in  the  Presence  of  Gradually  Changing 
Mass  Properties 

Table  7.7  shows  the  impact  of  the  aerodynamic  torques  on  the 
ability  of  the  mass  property  identification  algorithm  to  accurately 
estimate  gradually  changing  mass  properties.  During  the  identifi¬ 
cation  period  starting  at  t=70.36  sec,  the  filter  becomes  optimistic 
and  is  unable  to  converge  to  accurate  estimates. 


Table  7.7  New  Design  w/  Aerodynamics 
(Slow,  Continual  Changes  in  Mass  Properties) 

***  TIME  =  .52  SEC  *** 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX  IS: 

2248.20  -5.00  -63.60 

-5.00  1459.20  21.00 

-63.60  21.00  1178.60 

CENTER  OF  MASS  IS: 

83.75 

.12 

-.46 

UPDATED  MASS  PROPERTY  ESTIMATES: 

EST  INERTIA  MATRIX  IS: 

2243.05  -4.97  -63.47 

-4.97  1459.14  21.07 

-63.47  21.07  1178.50 

ESTIMATED  CENTER  OF  MASS  IS: 

83.67 

.04 

-.54 


MASS  PROPERTY  CHANGE  DETECTED  AT  T  =  41.32  SEC 


JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 


JET 

JET 

JET 

JET 

JET 

JET 

JET 

JET 


AT 

AT 

AT 

AT 

AT 

AT 

AT 

AT 


41.40  SEC 

41.60  SEC 

41.80  SEC 
42.00  SEC 
42.20  SEC 

42.40  SEC 

42.60  SEC 

42.80  SEC 
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Table 

7.7  (continued) 

JET 

#  6  AT 

T 

43.00 

SEC 

JET 

#  3  AT 

T 

= 

43.20 

SEC 

***  TIME 

=  43.28 

SEC 

*  *  ★ 

ACTUAL  MASS  PROPERTIES: 

INERTIA  MATRIX 

:  IS: 

2248.20 

-5. 

00 

•63.60 

-5.00 

1459. 

20 

21.00 

-63.60 

21. 

00  1178.60  j 

CENTER 

OF  MASS  IS: 

96.09 

12.46 

11.88 

UPDATED  MASS  PROPERT!^ 

•  ESTIMATES: 

EST  INERTIA  MATRIX 

IS: 

2233.49 

6. 

89 

■54.45 

6.89 

3457, 

09 

27.34 

-54.45 

27. 

34  1186.85 

ESTIMATED  CENTER  OF 

■  MASS 

IS: 

96.37 

10.97 

11.36 

MASS  PROPERTY  CHANGE 

DETECTED  AT  T  =»  7  0.36  SEC 

JET  FIRINGS  COMMANDED  BY  MASS  PROP  ID 

JET 

#  7  AT 

T 

- 

70.40 

SEC 

JET 

#  1  AT 

T 

= 

70.60 

SEC 

JET 

#  2  AT 

T 

= 

70.80 

SEC 

JET 

#  2  AT 

T 

= 

71.00 

SEC 

JET 

#  1  AT 

T 

= 

71.20 

SEC 

JET 

#  6  AT 

T 

= 

71.40 

SEC 

JET 

#  6  AT 

T 

= 

71.60 

SEC 

JET 

#  1  AT 

T 

= 

71.80 

SEC 

JET 

#  4  AT 

T 

= 

72.00 

SEC 

JET 

#  6  AT 

T 

72.20 

SEC 

JET 

#  6  AT 

T 

= 

72.40 

SEC 

JET 

#  4  AT 

T 

= 

72.60 

SEC 

JET 

#  1  AT 

T 

= 

72.80 

SEC 

JET 

#  4  AT 

T 

= 

73.00 

SEC 

JET 

#  6  AT 

T 

73.20 

SEC 

JET 

#  6  AT 

T 

- 

73.40 

SEC 

JET 

#  1  AT 

T 

- 

73.60 

SEC 

JET 

#  4  AT 

T 

= 

73.80 

SEC 

JET 

#  6  AT 

T 

31 

74.00 

SEC 

JET 

#  4  AT 

T 

3= 

74.20 

SEC 

JET 

#  6  AT 

T 

= 

74.40 

SEC 

JET 

#  1  AT 

T 

= 

74.60 

SEC 

JET 

#  4  AT 

T 

= 

74.80 

SEC 

(Firings 

continue 

until  end  of  simulation.  No  new  mass  property 

1  estimates  provided) 
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sliding  mode  control  law/linear  programming  jet  selection  algorithm 
combination  provides  somewhat  more  precise  tracking,  while  using 
less  fuel  (i.e.,  fewer  jet  firings),  than  the  existing  AFE  DAP  system 
consisting  of  a  proportional  control  law  and  a  table  look-up  jet 
selection  algorithm.  Also,  even  without  employing  the  mass  property 
identification  algorithm,  the  new  control  law/jet  selection  algorithm 
is  somewhat  more  robust  than  the  AFE  DAP,  Second,  the 
employment  of  the  mass  property  identification  algorithm  greatly 
extends  the  range  of  vehicle  mass  properties  over  which  the  attitude 
control  system  is  robust. 

The  results  of  the  simulation  tests  described  in  Section  7.4 
show  the  ability  of  the  candidate  control  system  to  not  only  make 
very  accurate  initial  mass  property  estimates,  but  to  also  detect  and 
adapt  to  mass  property  changes  that  occur  during  the  aerocapture 
maneuver  that  would  normally  cause  the  system  to  become  unstable. 
The  mass  property  identification  algorithm  can  detect  instantaneous 
or  gradual  mass  property  changes  by  comparing  the  measured  rate 
changes  due  to  the  jet  firings  to  predicted  rate  changes.  After  a 
change  is  detected,  the  operation  scheme  described  in  Chapter  6 
allows  the  control  system  to  continue  to  provide  attitude  control 
while  mass  property  identification  is  performed.  At  the  end  of  an 
identification  period,  the  updated  mass  property  estimates  are 
provided  to  both  the  control  law  and  the  jet  selection  algorithm 
enabling  the  control  system  to  continue  to  provide  stable,  efficient, 
attitude  control. 

Test  results  presented  in  Section  7.5  show  that  the  control 
system  is  robust  to  the  aerodynamic  torques  experienced  by  the 
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Chapter  8 

Conclusions  and  Recommendations 


8.1  Conclusions 

The  main  objective  of  this  research  has  been  the  design  of  a 
candidate  control  system  capable  of  providing  attitude  control  for  an 
aerocapture  vehicle  with  a  wide  range  of  possible  mass  properties. 
The  task  demanded  of  the  attitude  control  system  has  been  to  track 
the  desired  bank  angle  trajectory  commanded  by  the  guidance  law 
while  maintaining  trim  angle-of-attack  (a)  and  sideslip  angle  (3), 

Additionally,  the  control  system  was  required  to  be  able  to  perform 

this  task  over  a  wide  range  of  vehicle  mass  properties. 

The  candidate  control  system  consists  of  the  following  three 
main  elements:  a  nonlinear  sliding  mode  control  law,  an  adaptable 
linear  programming  jet  selection  algorithm,  and  a  mass  property 
identification  algorithm  based  on  a  second-order  extended  Kalman 
filter  design.  The  jet  selection  and  mass  property  identification 
algorithms  hav^e  both  been  used  previously  in  similar  works 

[20,21,24],  The  sliding  mode  control  law,  on  the  other  hand,  was 

designed  for  this  effort,  while  the  main  contribution  of  this  work  has 
been  the  integration  of  these  elements  into  the  candidate  design. 

The  integrated  candidate  contioi  system  design  does  meet  the 
main  objective  described  above.  Two  significant  results  were  de¬ 
monstrated  in  the  simulation  tests  discussed  in  Section  7.3.  First,  the 
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vehicle  during  the  aerocapture  maneuver.  As  the  spacecraft 
descends  deeper  into  the  atmosphere,  however,  the  aerodynamic 
torques  experienced  degrade  the  accuracy  of  the  mass  property 
identification.  At  lower  altitudes  where  the  magnitude  of  the 
aerodynamic  torques  are  greatest,  the  filter  is  not  able  to  converge  to 
accurate  mass  property  estimates. 

8.2  Recommendations  for  Future  Work 

Future  work  should  address  the  effects  of  aerodynamic  torques 
on  the  ability  of  the  mass  property  identification  algorithm  to 
accurately  estimate  mass  properties.  Additionally,  since  this  thesis 
addresses  only  a  subset  of  the  ultimate  goal  (described  in  Section 
1.5)  of  a  generic  attitude  control  system  for  use  with  any  aerocapture 
mission  and  any  vehicle,  future  research  should  concentrate  on 
achieving  the  ultimate  goal.  Specifically,  the  recommendations  are: 

(1)  Develop  an  estimation  algorithm  to  continually  identify 
the  aerodynamic  forces  and  torques  on  the  aerocapture 
vehicle.  With  this  knowledge,  the  effects  of  these  forces 
and  torques  can  be  subtracted  out  of  the  calculation  of 
the  residual  vector  so  that  the  mass  property 
identification  algorithm  can  more  accurately  estimate  the 
mass  properties  in  the  presence  of  significant  aero¬ 
dynamic  torques. 

(2)  Extend  the  mass  property  identification  algorithm  to 
accommodate  non-ideal  attitude  and  rate  measurements. 
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(3)  Test  the  candidate  system  in  a  simulation  with  an 

advanced  guidance  law  that  varies  both  bank  angle  and 

angle-of-attack. 

(4)  Develop  algorithms  to  identify  failed  jets,  along  with 

reduced  or  misaligned  jet  thrust. 

(5)  Develop  an  algorithm  to  identify  the  effectiveness  of 

alternate  types  of  control  actuators  (i.e.,  aerosurfaces  and 
control  moment  gyros)  so  that  the  attitude  control  system 
structure  can  be  applied  to  a  vehicle  employing  blended 
control. 

(6)  Apply  the  candidate  attitude  control  system  to  a  vehicle 
performing  a  Martian  aerocapture  to  ensure  the  cand¬ 
idate  system  is  viable  in  greatly  different  environments. 

Finally,  this  candidate  control  system  structure  is  not  limited  to 
attitude  control  of  aerocapture  vehicles.  It  may  be  applicable  to 
many  types  of  control  problems  for  which  parameters  impacting 
actuator  effectiveness  are  not  known  and  may  vary  significantly 
during  system  operation.  Future  work  should  investigate  the 
benefits  of  applying  this  control  system  structure  to  other  types  of 
control  problems. 
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Appendix  A 


Velocity  Angles  and  the  Body  Axes 

A.l  123  Axes,  Body  Axes,  and  Velocity  Angles 


Figure  A.l  123  Reference  System 


Figure  A.l  illustrates  the  orientation  of  a  local-vertical  local- 
horizontal  axes  system  denoted  as  123  with  respect  to  the  Earth 
fixed  axes  XYZ.  The  angle  of  longitude  the  angle  of  latitude  p,  and 
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the  radial  distance  iz  from  the  center  of  Earth  define  the  location  of 
the  center  of  the  123  axes  system  corresponding  to  the  center  of 
mass  of  the  aerocapture  vehicle.  The  1  axis  points  towards  the  north, 
the  2  axis  towards  the  east,  and  the  3  axis  towards  the  center  of 
Earth. 


Figure  A. 2.  Euler  Angles 


The  Euler  angles  0,  9,  and  \|r  shown  in  Figure  A. 2  define  the 
orientation  of  the  xyz  vehicle  body  axes  with  respect  to  the  123  axes. 
The  order  of  rotations  performed  to  achieve  the  body  axes  is: 


0 


\)/  about  the  3  axis 
0  about  the  2’  axis 
(p  about  the  1"  axis 

Bank  angle,  angle  of  attack,  and  sideslip  angle  (i.e.,  <(),  a,  and  P) 
are  considered  the  velocity  angles.  These  angles  determine  the 
orientation  of  the  body  axes  with  respect  to  the  velocity  vector. 
Figure  A.3  shows  how  the  velocity  angles  are  defined  with  respect  to 
the  body  referenced  axis  system. 


The  order  of  rotations  required  to  orient  the  body  axes  with  respect 
to  the  velocity  vector  is: 


<t)  about  V 

P  about  the  Zq  axis  (other  conventions 

use  P  about  -Zq) 

a  about  the  y  axis 


A.2  Rotation  Matrix  for  Conversion  From  Body  Rates 
to  Velocity  Angle  Rates 

Follov.'ing  the  examples  in  chapter  2  of  reference  [39],  the 
rotation  matrix  required  to  convert  body  rates  into  velocity  angle 
rates  can  be  calculated  with  two  rotations.  The  starting  point  is  an 
axis  system  ABC  with  the  A  axis  along  the  velocity  vector.  Since  the 
velocity  angles  a  and  P  are  independent  of  the  bank  angle,  the  first 
step  is  to  rotate  an  angle  p  about  the  C  axis. 
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The  ABC  axis  system  can  be  written  in  terms  of  the  new  x’y'z'  axis 
system  by 


A 

cos  P  -sin  P 

0 

x’ 

B 

= 

sin  p  cos  P 

0 

y' 

-C. 

0  0 

1 

z' 

(A.l) 


The  x'y'z'  axis  system  is  now  rotated  about  the  y’  axis  through  an 
angle  of  a  and  the  x'y’z'  can  be  written  in  terms  of  the  body 
referenced  axes  xyz. 


X 

y’ 


cos  a 
0 

L  -sin  a 


0 

1 

0 


sin  a 
0 

cos  a. 


(A.2) 


Combining  Equations  A.l  and  A. 2,  the  ABC  axis  system  can  be  written 
in  terms  of  the  body  referenced  system  from 
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A 

'  x  ' 

B 

=  R 

y 

-C- 

z 

(A.3) 


where  the  rotation  matrix,  R,  is  given  by 


R  = 


cos  a  cos  P 
cos  a  sin  ^ 
-sin  a 


-sin  ji 
cos  P 
0 


sin  a  cos  P 
sin  a  sin  p 
cos  a 


(A.3) 


The  resulting  rotation  matrix  can  be  used  to  express  velocity  angle 
rates  in  terms  of  body  angle  rates  and  vice  versa.  In  equation  form, 


’0' 

■p' 

ct 

=  R 

q 

r 

LpJ 

where  p,  q,  and  r  are  the  roll  rate  (about  the  x  axis),  pitch  rate  (about 
the  y  axis),  and  yaw  rate  (about  the  z  axis),  respectively. 
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Appendix  B 

Estimates  and  Bounds  for  Control  Law 
Design 


B.l  Inertia  Matrix 

Reference  [36]  provides  the  mass  properties  of  the  AFE  vehicle. 
At  atmosphere  entry,  the  inertia  matrix  is  predicted  to  be 


1 

■  2290.7 

-5.0 

-63.6  ■ 

1  = 

-5.0 

1490.9 

21.0 

slugs  • 

ft^ 

.  -63.6 

21.0 

1189.5  - 

atmosphere 

exit  (after  fuel 

consumption),  I 

is  predicted 

■  2196.2 

-5.0 

-63.6  ■ 

1  = 

-5.0 

1420.4 

21.0 

slugs  • 

ft^ 

.  -63.6 

21.0 

1165.2  J 

Assuming  these  estimates  are  accurate  to  within  +/-  10%,  the  range 
on  I  for  the  nominal  AFE  vehicle  is  assumed  to  be 


■  1976.6  -5.5  -70.0  ' 

'2519.8 

-4.5 

-57.2  ■ 

-5.5  1278.4  18.9 

<I< 

-4.5 

1640.0 

23.1 

.  -70.0  18.9  1048.7. 

.  -57.2 

23.1 

1308.5  - 
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Taking  the  midpoint  of  the  ranges  for  each  element,  the  estimated 
inertia  matrix  to  be  used  in  the  controller  design  is 


2248.2 

-5.0 

-63.6 

-5.0 

1459.2 

21.0 

-63.6 

21.0 

1178.6- 

slugs  •  ft^ 


(B.2) 


B.2  Input  Matrix 

In  order  to  calculate  the  estimate  and  bounds  for  the  input 
matrix  B  defined  by  Equation  3.30,  the  bounds  on  the  matrix  D 
(defined  in  Equation  3.35)  must  first  be  found.  Using  the  inertia 
matrix  estimates  and  bounds  found  above  and  assuming  the 
maximum  value  of  the  rotation  matrix  R  defined  by  Equation  A, 3  is 
the  3  by  3  identity  matrix,  from  Equation  3.45, 


-0.1080  -0.0004  -0.0048 
-0.0005  -0.1102  0.0001 
.-0.0096  0.0003  -0.0995- 


<  A< 


0.1380  0.0006  0.0076 
0.0008  0.1414  -0.0002 
.  0.0153  -0.0005  0.1244  . 


(B.3) 


With  B  equal  to  the  3  by  3  identity  matrix  as  specified  in  Equation 
3.45,  the  range  on  the  input  matrix  is  found  from  Equation  3.30, 


0.8920  -0.0004  -0.0048 
-0.0005  0.8898  0.0001 
.-0.0096  0.0003  0.9005  . 


<B< 


1.1380  0.0006  0.0076 
0.0008  1.1414  -0.0002 
-0.0153  -0.0005  1.1244  . 


(B.4) 
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B.3  Aerodynamic  Disturbances 

In  Equation  3.34,  the  system  of  equations  used  to  calculate  the 
control  law  gains  k,  dimax  is  the  maximum  angular  acceleration  in  the 
i*h  dimension  due  to  the  unmodelled  aerodynamic  torques.  Since  a 
and  P  errors  cause  torques  that  drive  the  vehicle  back  to  trim,  the 
main  aerodynamic  acceleration  of  concern  is  the  bank  acceleration 
due  to  velocity  angle  errors. 

Following  the  analysis  provided  in  reference  [37],  P  can  be 
assumed  small.  The  banking  moment  is  then  comprised  of  moments 
about  the  x  and  z  body  axes.  A  bank  moment  coefficient  can  be 
defined  as 


Q  =  Cn  sin  a  +  Cl  cos  a 


(B.5) 


where 

Cl'  =  bank  moment  coefficient 
Cn  =  yaw  moment  coefficient 
Cl  =  roll  moment  coefficient 


Assuming  a  and  P  remain  within  2  deg  of  their  trim  values,  from 
reference  [37] 


Cnmax  —  0.00562 
Clmax  =  0.00034 


The  maximum  bank  moment  coefficient  is  then  calculated  at  a  =  19 
deg  to  be 
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(B.7) 


C.  =  0.0005 

ifflax 

The  maximum  bank  moment  T())  is  calculated  from 


T4,  =  C’l  S  b  q 


(B.8) 


where 

S  =  reference  area  =  151.74  ft^ 
b  =  reference  length  =  13.99  ft 
q  =  dynamic  pressure  in  Ibs/ft^ 


Dynamic  pressure  is  a  function  of  velocity  and  atmospheric 
pressure. 


(B.9) 


where  p  is  the  atmospheric  density  and  v  is  the  magnitude  of  the 
velocity.  Figure  B.l  is  a  plot  of  q  vs.  altitude  during  the  simulation 
tests  performed.  Assuming  a  maximum  q  of  34.0  Ibs/ft^,  from 
Equation  B.8, 


=  36.08  Ibs  ft  (B.IO) 

The  disturbance  acceleration  due  to  this  aerodynamic  bank  moment 
is  estimated  to  be 


di 


max  — 


^XXr 


1976.6 


=  0.018  tads 


sec^ 


(B.ll) 
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The  maximum  disturbance  in  each  dimension  is  set  to  di. 


d:  =  0.018  sads 

‘max  o 

sec^ 


(B.12) 


altitude  (feet) 

Figure  B.l.  Dynamic  Pressure  vs.  Altitude 
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Appendix  C 

Simulated  Changes  to  Vehicle  Inertia 
Matrix 


This  appendix  describes  how  inertia  matrix  changes  due  to 
movements  of  the  center  of  mass  location  from  its  nominal  position 
were  calculated  during  simulation  testing. 

First,  center  of  mass  movements  were  simulated  without 
adding  mass  to  the  vehicle.  Assuming  a  point  mass  equal  to  half  the 
vehicle  mass  exists  at  the  nominal  center  of  mass  location,  the 
desired  displacement  of  the  center  of  mass  is  achieved  by  moving 
the  point  mass  to  twice  the  desired  displacement.  The  resulting 
center  of  mass  is  then  at  the  desired  location.  For  example,  for  the 
45  inch  center  of  mass  displacement  used  in  Section  6.4.2,  half  the 
vehicle  mass  was  moved  90  inches  from 


to 


cninominal 


83.60 

-0.03 

-0.61 


inches 


half  mass  position  = 


135.56 
51.93 
.  51.35  - 


inches 


(C.l) 


(C.2) 
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so  that  the  new  center  of  mass  location  is  45  inches  from  the  nominal 
location  at 


t*cninew 


The  new  inertia  matrix  elements  are  calculated  based  on  their 
definitions  given  in  reference  [30]; 

Ixx=i  mi(yf  +  z?)  (C.4) 

i=  1 

lyy  =  X  (Z?  +  Xf) 

i=  I 

kz  =  X  +  y?) 

i=  1 
N 

Ixy  =  -X  miXiyi 
i=  1 

N 

Ixz  ~  ■  ^  trii  X|  Z'l 
i=  1 

N 

Iyz  =  -X  miViZi 

i=  1 

where  the  vehicle  consists  of  N  mass  particles  each  with  mass  m;  and 
located  at  xiyizj.  Since  a  point  mass  was  removed  from  the  nominal 
center  of  mass  location,  the  inertia  due  to  the  remaining  mass  at  this 


109.58 
25.95 
L  25.37 


inches 


(C.3) 
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point  can  be  assumed  to  remain  equal  to  the  nominal  inertia  matrix. 
By  adding  half  the  mass  at  twice  the  desired  displacement  as 
described  above  (i.e.,  at  the  location  given  by  Equation  C.2  for  the  45 
inch  displacement  example),  the  new  inertia  matrix  elements  are 
calculated  by: 


Ixx  -  Ixxnominal  ^  (ym  +  ^m) 

lyy  =  lyynominal  ^  +  Xm) 

Izz  =  Izz„ommal  +  ^(xm  +  ym) 
Ixy  =  Ixynominal  "  ^  ym 
~  ^XZnominal  ^ 
lyz  ~  ^yZnominal  ” 


(C.5) 


where  m  is  the  total  vehicle  mass  and  Xm ,  ym .  and  Zm  are  the 
displacement  in  feet  of  the  half  mass  from  the  nominal  location  (e.g., 
Xm  =  [135.56  -  83.6J/12.0  =  4.33  ft  for  the  example  above). 

Using  Equations  C.5,  the  inertia  matrix  for  a  45  inch  center  of 
mass  displacement  is 


4554.31 

-1158.06 

-1216.66 


-1158.06 

3765.31 

-1132.06 


-1216.66 

-1132.06 

3484.71 


slugsft^ 


(C.6) 
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